特征波反演成像理論框架

王華忠,馮 波,王雄文,胡江濤,劉少勇,李 輝,周 陽

(波現(xiàn)象與反演成像研究組(WPI),同濟大學(xué)海洋與地球科學(xué)學(xué)院,上海200092)

特征波反演成像理論框架

王華忠,馮 波,王雄文,胡江濤,劉少勇,李 輝,周 陽

(波現(xiàn)象與反演成像研究組(WPI),同濟大學(xué)海洋與地球科學(xué)學(xué)院,上海200092)

地震波反演成像的核心問題是將解一個非線性(較)強的反問題轉(zhuǎn)化為提一個更凸的反問題并進行求解。存在強非線性性的主要原因是實測數(shù)據(jù)與要反演的模型參數(shù)之間的關(guān)系遠(yuǎn)非線性,其次是由于包含模型參數(shù)的控制方程不能很好地預(yù)測實際數(shù)據(jù)。因此,提出了特征波反演(characteristic wave inversion,CWI)成像的理論框架,基本思想是:不追求對實測波場全部波現(xiàn)象的模擬,而是模擬其中的部分特征波場;不一定追求波形逼近,但要盡可能利用到達時(相位)的逼近。相對于全波形反演(full waveform inversion,FWI),特征波反演由四個基本步驟組成:①特征波場(characteristic wave field,CWF)的提取;②波動理論的透射波層析成像;③波動理論的一次反射波層析成像;④最小二乘疊前深度偏移成像。特征波場提取是其中重要的環(huán)節(jié),包含三重含義:①波現(xiàn)象的分解(譬如,矢量波分解成標(biāo)量波以及一次波和多次波分解);②時空局部的、單震相的、帶方向的帶限波場的分解;③同相軸上地震子波的分解(譬如,提取子波的達到時、相位等)。特征波場提取基于壓縮感知的框架進行,其它三個線性化的參數(shù)反演環(huán)節(jié),首先考慮的是針對地下介質(zhì)參數(shù)層狀分布時的反射波反演成像,然后再考慮針對散射和繞射波的反演成像。數(shù)據(jù)域特征波反演在估計低波數(shù)速度信息時僅依賴同相軸上子波的到達時或/和相位信息,需盡量排除振幅對到達時和相位估計的影響。像域中的背景速度反演僅適宜基于到達時的反演,基于像的幅值的反演在理論上是不合理的。高波數(shù)參數(shù)估計時,首先進行方位角度反射系數(shù)的估計,在此基礎(chǔ)上進行散射強度的估計。CWI技術(shù)系列是推進經(jīng)典FWI走向?qū)嵱没恼_途徑,初步數(shù)值試驗結(jié)果證明了上述判斷。

反演成像;特征波場;非線性問題;凸問題;數(shù)據(jù)域反演;像域反演

巖性油氣藏的精確描述是勘探地震學(xué)的重要目標(biāo)。小尺度地質(zhì)體的精確構(gòu)造成像、振幅保真成像、反射系數(shù)甚至其它彈性參數(shù)的精確估計是勘探地震學(xué)的核心問題。隨著高性能計算機技術(shù)的快速進步以及全方位、高密度、寬頻帶的地震數(shù)據(jù)采集技術(shù)的飛速發(fā)展,以估計全或?qū)挷〝?shù)帶的彈性參數(shù)為目標(biāo)的全波形反演(FWI)方法再次成為勘探地震學(xué)的研究熱點。(各向異性)逆時偏移(RTM)方法的成功實用化則進一步推動了對FWI技術(shù)的研究。然而,由于疊前地震數(shù)據(jù)的不完備(如有限觀測孔徑、有限帶寬(尤其是缺失低頻)和空間不規(guī)則的觀測以及地表不一致的觀測,也包括強噪聲環(huán)境下的觀測)、正演模型的不完善(正演算子無法很好地預(yù)測實際波場中的各種波現(xiàn)象)、初始模型距離反演的“真”參數(shù)模型太遠(yuǎn)、地震子波的未知且空變等多種因素制約,使得經(jīng)典FWI難以達到其理論預(yù)計的目標(biāo),甚至在實際應(yīng)用中根本不收斂。

本質(zhì)上,經(jīng)典FWI是一個非線性性很強的反演問題,根本原因在于實際觀測的地震波場的變化與要反演參數(shù)的變化之間關(guān)系很弱。另外,包含要反演參數(shù)的控制方程不能很好地預(yù)測實測數(shù)據(jù)也是一個重要因素。而地震波反演成像的核心問題是將解一個非線性(較)強的反問題轉(zhuǎn)化為提一個更凸的反問題并進行求解。如何將一個非線性強的反問題轉(zhuǎn)化為凸性更好的反問題并進行求解,一直是地球物理反問題研究人員關(guān)注的焦點。典型地,包括引入波場的Laplace變換[1]、引入波場的對數(shù)變換[2]、引入波場的Hilbert變換[3-4],還有不進行數(shù)據(jù)逼近而利用像域中成像道集的拉平或聚焦建立目標(biāo)泛函的反演方法[5]。更為全面地討論如何構(gòu)建較凸反問題的,應(yīng)該是波場重構(gòu)反演(wavefield reconstruction inversion,WRI)方法[6-7],其基本思想是通過引入增廣的狀態(tài)變量解決正問題不能很好地預(yù)測實際地震數(shù)據(jù)中波現(xiàn)象的問題。但是,引入增廣的狀態(tài)變量并沒有明確地指出僅僅特征波場和特征參數(shù)之間才存在更為密切、更為線性的聯(lián)系。而勘探地震中地下介質(zhì)參數(shù)的分布和疊前數(shù)據(jù)本身都是具有特征的,特征波場和特征參數(shù)之間才存在更為密切的關(guān)系[8-9]。

本文提出僅僅用與特征參數(shù)關(guān)系密切的特征波場進行反演的思想與方法。特征波場與反演參數(shù)之間的關(guān)系更緊密,基于特征波場的反問題凸性更好。特征波反演成像由四個基本步驟組成:特征波場的提取、波動理論的透射波層析成像、波動理論的一次反射波層析成像、最小二乘疊前深度偏移成像。特征波反演成像的中心工作在于特征波場的提取。特征波反演成像中梯度的預(yù)條件要相對簡單,因為此時特征波對應(yīng)的梯度干擾大大減弱。另外,反射子波特征對于估計反射系數(shù)或波阻抗一類的反演成像十分重要[10]。特征波成像是基于對實測炮集中所包含的波現(xiàn)象認(rèn)識的一種反演框架。我們認(rèn)為,只有經(jīng)過相當(dāng)長時間的特征波反演成像方法研究,經(jīng)典FWI中所蘊含的問題才能逐步得到認(rèn)識和解決,進而發(fā)展出工業(yè)界實用化的反演成像技術(shù)。

1 求解強非線性反問題或提出更凸的反問題

地震波反演成像是勘探地震的核心問題,它是一個標(biāo)準(zhǔn)的非線性反演問題,其非線性程度隨著地下介質(zhì)的復(fù)雜程度提高(譬如小尺度速度異常體存在、橫向變速劇烈等)而提高。實測波場變化與地下介質(zhì)變化之間的非線性程度決定了地震波反演問題的非線性程度。針對勘探地震中的參數(shù)反演問題,我們是專注于尋找解決較強非線性反演問題的方法,還是專注于提一個凸性更強的反演問題來求解？

ROCKAFELLAR[11]指出,優(yōu)化問題的核心不在于問題的線性與非線性,而是問題的凸性與非凸性。顯然,提一個更凸的反問題比奮力解一個強非線性反問題更合理。

自TARANTOLA[12]建立以波動方程為基礎(chǔ)的地震反演理論近30年(尤其是最近幾年)來,隨著以“兩寬一高”為代表的地震數(shù)據(jù)采集技術(shù)以及以GPU為代表的高性能計算機技術(shù)的進步,FWI方法和技術(shù)取得了長足的進步。專家學(xué)者在降低反演問題非線性性方面做了很多有益的探索,提出了很多有效的方法?；舅枷胧菑臄?shù)據(jù)或波場入手來改進反演成像方法,典型的做法包括:從低頻數(shù)據(jù)到高頻數(shù)據(jù)的多尺度反演[13-14];從近偏移距數(shù)據(jù)到遠(yuǎn)偏移距數(shù)據(jù)的反演[15];從小到大調(diào)節(jié)Laplace阻尼參數(shù)進行的由初至?xí)r刻不斷放大時窗的反演[1];定義不同誤差泛函的反演[3,16-22]。另一類方法則是將傳統(tǒng)的數(shù)據(jù)域反演擴展到成像域的反演[5,23-24]。

上述反演方法并沒有從更為抽象的角度主動地將地震波反演問題提成一個更凸的問題。本文希望在數(shù)據(jù)的特征表達和參數(shù)的特征表達概念下,將經(jīng)典FWI轉(zhuǎn)化為一個臨近的、更凸的反問題進行求解。

2 Bayes反演理論的缺陷

(1)

式中:d代表觀測數(shù)據(jù),m為系統(tǒng)(或模型)參數(shù)。公式(1)表示參數(shù)估計的結(jié)果是由觀測數(shù)據(jù)所對應(yīng)的模型參數(shù)出現(xiàn)的概率所決定的。

由于后驗概率密度P(m/d)很難估計,所以將其轉(zhuǎn)化為較容易估計的模型空間的先驗概率密度函數(shù)P(m)、數(shù)據(jù)空間的概率密度函數(shù)P(d)和條件概率密度函數(shù)P(d/m)來表示,其中P(d/m)體現(xiàn)了正算子的不確定性。計算后驗概率密度函數(shù)P(m/d)的Bayes公式為:

(2)

(3)

式中:CD為數(shù)據(jù)協(xié)方差矩陣;CM為模型協(xié)方差矩陣;dobs代表實測數(shù)據(jù);dcal=g(m)代表正問題;mprior代表參數(shù)的先驗值。

(3)式存在的問題是隱含了一個很基本的假設(shè):正問題g(m)是可以解釋實測數(shù)據(jù)dobs的,不能解釋的殘差是符合高斯或某種廣義高斯分布的。這種假設(shè)與勘探地震的實際情況很不相符,必須進行修改。

TARANTOLA[25]指出:正問題不能很好地解釋實測地震數(shù)據(jù),其引入的數(shù)據(jù)殘差與觀測噪聲(非信號部分)引入的數(shù)據(jù)殘差可能是同一水平的。

我們認(rèn)為,考慮到正演模型選擇過于簡單(譬如僅考慮單參數(shù))及震源函數(shù)未知等因素的影響,正問題的表達不合適引起的數(shù)據(jù)逼近誤差很可能超過觀測噪聲引起的數(shù)據(jù)逼近誤差。這就使經(jīng)典的FWI完全失去了意義。Bayes估計理論似乎蘊含了正問題不準(zhǔn)確引起的問題,但要么忽視它,要么沒有具體的解決方案。

以下我們引入正問題的計算結(jié)果——狀態(tài)變量,并將狀態(tài)變量作為與模型參數(shù)同等地位的待反演參數(shù)看待。這是深化對地震波反演問題的認(rèn)識、將非線性的地震波反演成像問題提成一個更凸的反演問題的關(guān)鍵。

3 對狀態(tài)變量獨立的地震波反演問題的認(rèn)識

狀態(tài)變量是最優(yōu)控制理論[26]描述系統(tǒng)狀態(tài)的變量,與系統(tǒng)參數(shù)具有同等地位,共同滿足控制方程。

勘探地震中,引入狀態(tài)變量和系統(tǒng)參數(shù)的正問題的描述方法如下:將激勵、探區(qū)介質(zhì)分布、響應(yīng)視為一個系統(tǒng),對應(yīng)該系統(tǒng)的非線性正問題的數(shù)學(xué)描述可以抽象為:

A(m)u=-δ(x-xs)f(x,t)=-f(xs,t)

(4)

式中:A為正演算子;u為要正演合成的波場,即狀態(tài)變量(系統(tǒng)響應(yīng));f為震源函數(shù)(激勵函數(shù))。

在地震波反演中,A的一般形式為:

(5)

A的形式受所選擇的波動方程的影響,如常密度聲波/變密度聲波方程、各種擬聲波方程、彈性波方程等。

野外觀測到的地震數(shù)據(jù)可以表述為:

(6)

此處,P(·)表示針對實測數(shù)據(jù)和狀態(tài)變量的變換算子;η表示隨機噪聲函數(shù),在大部分情況下假設(shè)它是高斯白噪的。

由(6)式定義實測數(shù)據(jù)和狀態(tài)變量之間的關(guān)系,蘊含了引入變換算子使得實測數(shù)據(jù)和狀態(tài)變量進行匹配的思想,也可以說引入了讓狀態(tài)變量獨立出來作為一個待處理的參量看待的思想,不再認(rèn)為狀態(tài)變量完全由系統(tǒng)參數(shù)、控制方程和激勵函數(shù)決定(理論上的確由它們確定)。這樣的觀點對后續(xù)的反演成像非常有意義。

地震數(shù)據(jù)反演成像主要包括相互有緊密聯(lián)系的兩條途徑:一是Bayes估計理論;二是約束最優(yōu)化理論。廣義高斯分布假設(shè)下代價函數(shù)最小反問題優(yōu)化求解方法與約束最優(yōu)化理論基本上是等價的。下面我們從約束優(yōu)化理論展開狀態(tài)變量獨立的地震波反演成像問題的討論。

最簡單的FWI成像目標(biāo)函數(shù)定義為:

(7)

式中:dobs為野外觀測數(shù)據(jù),dcal(m)為正演模擬數(shù)據(jù);M是系統(tǒng)(模型)參數(shù)空間。(7)式就是數(shù)據(jù)殘差為高斯分布假設(shè)下的最大似然估計方法。另一個蘊含的假設(shè)是正問題的解能夠解釋實測地震數(shù)據(jù),并使二者的殘差滿足高斯分布。這與地震波反演成像的實際情況相差甚遠(yuǎn),反演結(jié)果有很強的多解性,完全不能滿足實際應(yīng)用的要求。因此,人們對(7)式進行了各種各樣的修改,發(fā)展了很多經(jīng)典的FWI變種形式。一般地,是沿著Bayes估計的路線提出各種修改方法。概念上,各種典型的修改方法可以抽象為:

(8)

4 狀態(tài)變量獨立的地震波反演問題的提法

狀態(tài)變量獨立的地震波反演問題可以提成如下控制方程約束下的優(yōu)化問題,即:

(9)

式中:Pcal(·)為施加在正演波場上的變換函數(shù);Pobs(·)為施加在觀測數(shù)據(jù)上的變換算子,它可以是非線性算子,一般定義為線性算子;U是狀態(tài)變量空間;A(m)u=q為控制(或狀態(tài))方程,q為激勵(震源)函數(shù)。狀態(tài)變量空間U和模型參數(shù)空間M構(gòu)成了滿足控制方程A(m)u=q的聯(lián)合空間,地震波反演成像就是在此聯(lián)合空間中尋找最優(yōu)解。事實上,激勵函數(shù)q∈Q也是需要反演的,但是在勘探地震中,一般假定激勵函數(shù)是已知的。關(guān)于模型參數(shù)和狀態(tài)變量的罰函數(shù)(約束函數(shù)),為了討論簡單計,沒有包括在上式中,但這不會影響下述問題的討論。

5 特征波反演方法

我們認(rèn)為,無論如何選擇控制方程、如何描述相應(yīng)的定解問題,都不可能完全解釋實測地震數(shù)據(jù)。既然如此,我們從實測數(shù)據(jù)中篩選出能被控制方程較準(zhǔn)確預(yù)測的波場,將這部分波場稱為特征波場(CWF),基于CWF發(fā)展一套特征波反演理論。特征地震波反演成像問題的凸性有明顯的改善,為簡潔計,下面的討論沒有考慮正則化問題。

基于CWF的特征波反演成像方法技術(shù)的核心在于特征波場數(shù)據(jù)的提取。實質(zhì)上,特征波場也對應(yīng)模型參數(shù)的特征成分。

壓縮感知框架下的特征波場提取問題可以提為:

(10)

在提取特征波場的基礎(chǔ)上,特征波場反演問題提為:

(11)

式中:R是采樣算子,AC為描述特征波傳播的波場正傳算子,uC為正演計算的特征波場,mC為與特征波場對應(yīng)的特征參數(shù)。它們共同滿足控制方程AC(mC)uC=q。由于我們僅僅考慮特征波場,因此,可以比較完整地進行特征波場的預(yù)測,(11)式可以重寫為:

(12)

VAN LEEUWEN等[7]提出如下波動方程約束的優(yōu)化反演問題(波場重構(gòu)反演方法):

(13)

上式也可以寫為:

(14)

可見波場重構(gòu)反演方法的思想也是提取實測波場中滿足控制方程的部分參與反演,即認(rèn)為u=P-1dobs,并利用(13)式定義了如下的增廣狀態(tài)變量:

(15)

但是,VAN LEEUWEN等沒有提出可行的計算重構(gòu)波場(或計算增廣狀態(tài)變量)的方法。

我們堅持用特征波場進行地震波反演的理念受到如下線性問題求解方式的啟發(fā):

(16)

特征波反演的核心問題之一在于解問題(10)式。相比于波場重構(gòu)反演方法[6],特征波場反演的物理意義十分明確,可操作性很強,據(jù)此可以構(gòu)建出有效的反演技術(shù)路線。尤其是,公式(10)的求解一般在壓縮感知框架下進行,可以將先進的信號處理理論運用于特征波提取及反演成像中。更重要的是,特征波的提取不僅僅是現(xiàn)代信號框架(或壓縮感知框架)下的操作,而是按照波動理論的物理意義進行操作(譬如僅僅提出透射波場和一次反射波場),盡可能建立起特征波場與特征參數(shù)之間的聯(lián)系,使得提出的特征波反演成像問題更凸。

特征波反演成像的基本過程可以描述如下。首先,我們將一般的地震數(shù)據(jù)表達成特征數(shù)據(jù),提取特征數(shù)據(jù)的過程可用下式抽象地表示:

(17)

在特征波場的概念下,一般的正問題可寫成如下積分形式:

δvb(x);xs,xr,x,t]mC(x)dx

(18)

式中:uC(xs,ps,xr,pr,t)是特征波場數(shù)據(jù),K代表靈敏度核函數(shù),mC(x)是特征模型參數(shù)(即反射界面的反射系數(shù));vb(x)是背景速度。在K的性質(zhì)比較好(譬如其逆算子存在)的情況下,特征波場數(shù)據(jù)與特征模型參數(shù)之間建立起了一對一的映射關(guān)系,此時的反演成像問題凸性一定是好的。

在上述正問題的基礎(chǔ)上,可以建立如下地震波反演成像問題:

(19)

式中:uC[vb(x)+δvb(x),mC(x)]為計算的特征波場,δvb(x)是背景速度的低波數(shù)擾動;mC(x)則為速度場的高波數(shù)擾動(反射系數(shù)),它代表了速度模型的幾何特征(即反射界面的幾何形態(tài)信息)。假定背景速度vb(x)是已知的,δvb(x)和mC(x)是要反演的模型參數(shù)。事實上,背景速度vb(x)也要通過透射波層析反演來確定。

成像域的反射波層析也可以基于特征波反演框架進行描述。事實上,在特征波場情形下,成像域的模型參數(shù)層析反演應(yīng)該與數(shù)據(jù)域的模型參數(shù)層析反演等價。

如果進行旅行時層析,無論是數(shù)據(jù)域還是成像域,基本的問題都是時差測量問題。數(shù)據(jù)域中,實測數(shù)據(jù)同相軸到達時是確定的,實測數(shù)據(jù)和合成數(shù)據(jù)對應(yīng)同相軸的時差測量結(jié)果可以認(rèn)為是絕對的;成像域中,反射界面真深度是未知的,基于成像道集的剩余深度差(residual moveout,RMO)測量是相對的。這也是為什么需要在時間域進行層析速度反演的重要原因。實測數(shù)據(jù)和合成數(shù)據(jù)對應(yīng)同相軸的時差測量是當(dāng)前層析速度反演的一個核心問題,解決了此問題,FWI中的周期跳躍問題就不存在了。AWI(Adaptive Waveform Inversion)方法[22,29]的本質(zhì)也是試圖解決這個問題。

6 特征波場提取的思想與方法

從實測地震數(shù)據(jù)中根據(jù)反演成像的要求提取出特征波場的問題可以抽象為如下變換:

(20)

其中,變換算子P(·)可以根據(jù)反演成像方法的需要選擇各種形式,關(guān)鍵是構(gòu)建特定的反演方法時想利用什么樣的波現(xiàn)象、想利用該波現(xiàn)象的什么特征。譬如,對波場的Laplace變換[1]、對波場的Hilbert變換[3]、對波場的對數(shù)變換[27,30]等。這些變換算子都側(cè)重于提取反射子波的特征。

我們提取的特征波場首先考慮的是區(qū)別利用不同的波現(xiàn)象,然后再進一步考慮利用特定波現(xiàn)象的更進一步的特征。我們認(rèn)為,特定波現(xiàn)象的方向傳播特征是很基本的波場特征,因此,特征波場變換算子至少包含三重含義:①區(qū)分波現(xiàn)象;②提取時空局部的方向波場;③更進一步的子波波形特征的提取。

區(qū)分波現(xiàn)象主要是分開縱橫波、透射波和反射波,繞射波和反射波,還包括分開一次波和多次波。提取時空局部的方向波場主要有基于基函數(shù)族和基于數(shù)據(jù)驅(qū)動的兩種線性局部平面波提取方法。子波特征的提取主要是提取帶限子波的到達時和相位。

這三個方面目前都存在問題,尤其是在實測地震數(shù)據(jù)有較強噪聲的情況下。帶限子波的到達時和相位因為與振幅攪在一起,其估計從概念上就存在問題。因此,帶限子波波形特征的表達和提取對于反演成像而言是更為重要的研究內(nèi)容,否則可能連到達時的含義都說不清楚,速度層析反演的精度完全無法評價。子波振幅的影響因素就更復(fù)雜了。

本文主要關(guān)注的是在壓縮感知框架下如何提取時空局部的方向波場。這個問題被提為如下變換域的反問題:

(21a)

式中:A是變換矩陣,列向量構(gòu)成基函數(shù)族。數(shù)據(jù)的特征表達原則上也可以用如下數(shù)據(jù)驅(qū)動的方式:

(21b)

式中:‖X‖*代表核范數(shù),是矩陣XSVD分解后的奇異值的和;ARφH為數(shù)據(jù)變換算子,代表將SVD分解后的奇異值恢復(fù)成與實測數(shù)據(jù)采樣一致的數(shù)據(jù),其中R是采樣算子,φH是將SVD分解后的奇異值恢復(fù)到原數(shù)據(jù)域的算子。數(shù)據(jù)驅(qū)動的地震數(shù)據(jù)特征表達是目前地震數(shù)據(jù)預(yù)處理的一個研究熱點[31-33]。相對于基于基函數(shù)的波場特征表達,數(shù)據(jù)驅(qū)動的方法物理意義更難把握,因為特征向量的物理意義不是很明晰。

圖1展示了抽取特征波場的基本思想和做法。其中圖1a為匹配追蹤方法抽取特征波場的示意圖，從左到右包含4列：第1列為包含3個同相軸(①～③)的局部波場；第2列為匹配追蹤方法抽取水平同相軸②后的結(jié)果；第3列為再抽取斜同相軸③后的結(jié)果；第4列為3個同相軸全部抽取后的結(jié)果。圖1b為圖1a對應(yīng)的τ-p譜。圖1c展示了數(shù)據(jù)恢復(fù)情況：第1列為原數(shù)據(jù)(圖1a第1列所示)與反變換數(shù)據(jù)的差；第2列展示了恢復(fù)出的水平同相軸②；第3列展示了恢復(fù)的斜同相軸③；第4列展示了恢復(fù)的斜同相軸①。至此，比較完整地恢復(fù)出了原數(shù)據(jù)，實現(xiàn)了抽取和恢復(fù)的全過程。發(fā)展該方法的另一個目的是解決層析反演中實測數(shù)據(jù)同相軸和合成數(shù)據(jù)同相軸對應(yīng)的時差測量問題，試圖通過這樣的方法避開周期跳躍現(xiàn)象。

圖1 特征波場的抽取a 匹配追蹤方法抽取特征波場方法的示意; b 與抽取的特征波場相對應(yīng)的τ-p譜; c 數(shù)據(jù)恢復(fù)

7 數(shù)據(jù)域的特征波反演方法

數(shù)據(jù)域地震波反演的經(jīng)典方法自然是FWI方法。經(jīng)典FWI的目標(biāo)是估計全或?qū)挷〝?shù)帶的參數(shù)場。但是,到目前為止的實踐表明,即使用經(jīng)典FWI結(jié)合低頻長偏移距數(shù)據(jù)估計低波數(shù)(背景)的速度，也不是任何情形都能收斂的。這正是經(jīng)典FWI衍生出很多變種的原因,但這也造成了FWI概念的濫用。

此處我們要討論的是基于CWF的特征波場反演(CWI),包括四個關(guān)鍵步驟:①用壓縮感知方法提取所要的特征波場(目前主要是初至波、早至波、一次反射波,提取特征波場的目的主要是用來進行層析成像,估計背景參數(shù)的變化);②基于到達時(相位)的透射波波動理論層析反演;③基于到達時(相位)的一次反射波波動理論層析反演;④主要針對方位角度反射系數(shù)的最小二乘疊前深度偏移成像。

特征波場必須與要反演的參數(shù)之間有更緊密的關(guān)系才能使得構(gòu)建的反演問題的凸性更好。我們提出如下的模型參數(shù)分解方法:

(22)

式中:R(x,θ,φ)代表方位角度反射系數(shù),θ代表入射角度,φ代表方位角。

基于上述模型分解,波場也可以分解為:

(23)

等式右端三項的具體含義分別為:第一項對應(yīng)透射波;第二項對應(yīng)一次反射波;第三項對應(yīng)二次散射波。更高階的散射波沒有列出。地下存在散射點引起的場的擾動可以表示為:

(24)

式中:δG1(xr|xs)代表一階擾動場;δG2(xr|xs)代表二階擾動場;更高階的擾動場被忽略。當(dāng)弱散射假設(shè)成立時,更高階的擾動場很弱,忽略掉是合理的。

對應(yīng)地,零階擾動波場(透射波場)、一階擾動波場和二階擾動波場的波場擾動相對于介質(zhì)擾動的梯度分別定義如下。

(25a)

xs)δm(ξ)dξ=-ω2G0(xr|x)G0(x|xs)

(25b)

G0(η|ξ)δm(ξ)G0(ξ|xs)dξdη

δm(ξ)G0(ξ|x)G0(x|xr)dξ

=-ω2G0(x|xs)δG1(x|xr)+G0(x|xr)ω4·

=-ω2G0(x|xs)δG1(x|xr)-ω2G0(x|xr)·

δG1(x|xs)

(25c)

(25a)至(25c)式分別定義了聯(lián)系波場擾動和介質(zhì)參數(shù)擾動的積分核(Kernel)。波動理論下,3個梯度項的示意圖分別如圖2a,圖2b和圖2c所示。

更進一步地,對于一階擾動和二階擾動而言,必須充分利用地下介質(zhì)的既有特征。我們將油氣地震勘探面對的地下介質(zhì)的特征表述為:在空間上廣泛分布的層狀沉積層,加上火山活動、構(gòu)造運動等導(dǎo)致的大小尺度速度異常體(或波阻抗變化的異常體)。即地下介質(zhì)的分層性是勘探地震所面對的介質(zhì)的重要特征,產(chǎn)生繞射和散射的介質(zhì)特征是次要的。

因此,二階擾動波場的層析成像和一階擾動波場的偏移成像都利用了這樣的特征。只有對反射界面而言,才會產(chǎn)生圖2c所示的“兔子耳朵”響應(yīng)。對于層析成像而言,可以利用反偏移計算圖2c所示的二階擾動響應(yīng);對于偏移成像而言,可以針對方位角度反射系數(shù)進行成像。

實質(zhì)上,FWI是一種假設(shè)地下介質(zhì)為背景+散射體(點)模式的反演成像方法。這樣的成像方法應(yīng)該對應(yīng)封閉面上全方位激發(fā)接收的觀測方式,而不是勘探地震僅在地表面進行的觀測方式。我們認(rèn)為這也是FWI方法或逆散射成像方法不能很好地用于油氣地震勘探領(lǐng)域的重要原因。

基于上述分析,我們認(rèn)為,將FWI分為如CWI所闡述的三個線性反演步驟在數(shù)學(xué)物理上是合理的。更進一步地,我們認(rèn)為背景速度和各向異性參數(shù)與一階擾動和二階擾動場中的同相軸的到達時(相位)關(guān)系密切,這些參數(shù)的估計不僅要利用到達時(相位),還要盡可能剝離振幅對達到時(相位)測量的影響。事實上,層析反演的理論框架基本是定型的,實際數(shù)據(jù)到達時(相位)、合成數(shù)據(jù)到達時(相位)以及對應(yīng)同相軸上二者之間的差的測量是層析成像方法與技術(shù)的核心問題。我們提出壓縮感知框架下的CWF提取,其中一個重要的目的就是發(fā)展出有效的對應(yīng)同相軸上實際數(shù)據(jù)到達時(相位)和合成數(shù)據(jù)到達時(相位)之間的差的測量方法。

圖2 透射波擾動響應(yīng)(a)、一階擾動場響應(yīng)(b)和二階擾動場響應(yīng)(c)

這是我們堅持的將FWI實用化的具體路線。鑒于篇幅,具體數(shù)值結(jié)果見參考文獻[10]。

8 成像域的特征波反演

地震波反演也可以在成像域中進行,基本思想是不同角度(或偏移距)的數(shù)據(jù)產(chǎn)生的像理論上應(yīng)出現(xiàn)在相同的空間位置上,即所謂的正確的成像速度產(chǎn)生正確的像,錯誤的成像速度產(chǎn)生錯誤的像。也可以說成像速度的擾動產(chǎn)生了像的擾動,成像速度的擾動與像的擾動之間存在線性關(guān)系,建立這樣的線性關(guān)系就可以進行成像速度的反演。正確的像是對同一地下像點而言,不同角度或偏移距的照明數(shù)據(jù)產(chǎn)生的像出現(xiàn)在該地下成像點上,表現(xiàn)為成像道集拉平或聚焦在一起。我們認(rèn)為,這與Bayesian框架下的數(shù)據(jù)域的參數(shù)估計理論有著本質(zhì)區(qū)別,具體表現(xiàn)在:①在成像域僅僅能估計光滑的成像速度,因為高波數(shù)的速度變化對產(chǎn)生拉平或聚焦的道集并無貢獻,而道集的拉平或聚焦正是建立像域反演目標(biāo)泛函的基礎(chǔ);②像域的成像速度估計缺乏概率統(tǒng)計意義下的解釋。

特征數(shù)據(jù)體產(chǎn)生的像可以用下式表示:

(26)

構(gòu)建像域的反問題有兩條途徑:一是通過檢測成像道集的拉平或聚焦情況建立目標(biāo)泛函,通過對目標(biāo)泛函的求導(dǎo),建立梯度導(dǎo)引類的優(yōu)化反演方法,估計成像速度[5,34-37];二是直接建立像的擾動與成像速度擾動之間的線性關(guān)系[38-39],基于此線性關(guān)系進行成像速度的層析反演。

首先討論第一種方法。SHEN[35]給出了兩類所謂的差異相似性優(yōu)化(differential semblance optimization,DSO)目標(biāo)函數(shù),分別描述地下偏移距(h)道集在零偏移距處的聚焦特性和角度(θ)道集的拉平特性。目標(biāo)函數(shù)分別為:

式中:P(h)=h為零階差分算子;P(θ)=?/?θ為一階差分算子。

對于公式(27),當(dāng)速度正確時,地下偏移距道集聚焦于h=0處,目標(biāo)函數(shù)取得極小值。該目標(biāo)函數(shù)的凸性較好,且不受周期跳躍問題制約,但對中小尺度的速度擾動不敏感,速度估計的精度低。SHAN等[40-41]指出,(27)式的收斂受到成像道集振幅的影響。當(dāng)?shù)叵抡彰髂芰枯^弱時,需要引入照明補償?shù)?否則強反射層的像會主導(dǎo)泛函能量,影響泛函收斂。此外,(27)式的梯度中含有一些近垂直的條帶狀假象,影響反演的收斂。FEI等[42]通過修改成像差(用關(guān)于h的一階差分算子作用于地下偏移距域的像)的定義,改善泛函梯度的性質(zhì)并加速反演收斂。

根據(jù)成像原理,當(dāng)速度正確時,角度道集的像是拉平的,相鄰角度像的差異較小。由此可以得出公式(28)描述的角度道集DSO目標(biāo)函數(shù)。由于(28)式只需要對相鄰角度的像進行差分運算,所以弱化了周期跳躍問題的制約,降低了對初始模型的精度要求。然而,差分運算會放大角度道集中高頻變化的影響,對角度道集的質(zhì)量也有一定的要求。

SOUBARAS等[28]提出了角度道集疊加能量最大化方法,目標(biāo)函數(shù)可以寫為:

(29)

式中:I(x,z,θ)為角度域成像道集。該方法的精度高于上述DSO類的反演方法。然而,當(dāng)初始速度與真實速度的差異較大時,角度道集的RMO較大,甚至大于半個子波長度。此時,若按照公式(29)對所有角度疊加,會遇到與FWI類似的周期跳躍問題。因此,該方法對初始模型的要求較高。

值得提及的是,地下偏移距成像道集是不區(qū)分反射和散射(繞射)的,都可以產(chǎn)生對應(yīng)當(dāng)前偏移速度的聚焦道集;角度成像道集是針對反射波的,只有反射界面才有對應(yīng)的角度成像道集。因此,(27)式定義的泛函和(28)及(29)式定義的泛函適用于不同的介質(zhì)參數(shù)分布情形。我們認(rèn)為,反射界面的存在是勘探地震的特征,應(yīng)該充分利用這樣的特征。

在上述目標(biāo)泛函特性分析的基礎(chǔ)上,可以進入梯度導(dǎo)引的反演過程中。但是,我們認(rèn)為利用地下偏移距聚焦成像道集,基于聚焦好壞建立目標(biāo)泛函進行成像速度的反演理論上是存在問題的。不能用聚焦好壞(能量聚焦度)作為目標(biāo)函數(shù)進行優(yōu)化反演、更新成像速度。原因之一是目標(biāo)泛函的梯度不能指示速度修正量的正或負(fù)(這與數(shù)據(jù)逼近泛函不同);原因之二是影響振幅的因素太多,聚焦能量最小并不一定預(yù)示背景速度的修正是正確的。角度道集疊加能量最大化目標(biāo)泛函同樣存在不合適地利用振幅的問題。利用角度道集之間的像差建立目標(biāo)泛函進行速度反演是合理的邏輯,但是其中AVA關(guān)系的影響及像的振幅畸變的影響會導(dǎo)致方法不穩(wěn)定?？傊?我們不認(rèn)為在像域進行上述基于像的振幅構(gòu)建泛函、進行背景速度反演的方式是正確的選擇。DSO之類的像域背景速度估計概念20多年來始終無法在勘探地震中實際應(yīng)用很能說明這一問題。

以下討論構(gòu)建像域反問題的第二種方法。首先建立像的擾動與成像速度擾動之間的關(guān)系。由(26)式可見,像的擾動必然是其中波場的擾動引起的,而場的擾動又是由于成像速度的擾動(即成像速度與成像真速度之間的差)引起的。因此,引入像的擾動為:

(30)

很顯然,場的擾動與速度擾動之間的關(guān)系是由散射場的Born計算公式?jīng)Q定的,即:

(31a)

(31b)

將(31)式代入(30)式可以導(dǎo)出像的擾動和成像速度擾動之間的線性關(guān)系:

(32)

用共軛梯度法求解上式可以得到更新的成像速度。公式(30)至公式(32)概念上是沒有問題的,抽象地講,這些公式可以用于不區(qū)分散射和反射的一般情形,但是具體應(yīng)用時要區(qū)分是用于反射情形還是用于散射情形。在實際應(yīng)用中,(32)式右端項表示的像的擾動很難有效獲取,因為真像是未知的！在針對反射情形的角度成像道集中,像的擾動代表道集沒有拉平,本次偏移速度得到的像與真像之間有差異,并且最主要的差異就是成像深度差RMO！因此,還是要利用特征波場(局部平面波場)與特征像(角度成像道集)的概念,建立起它們之間的聯(lián)系,進行基于RMO的波動理論的角度成像道集層析成像。其中,有效的RMO的檢測是最為核心的問題。在此基礎(chǔ)上,有必要和有可能的話,再基于公式(30)至公式(32) 進行層析成像。這是我們一直堅持的像域中速度估計與建模的基本理念。鑒于篇幅,具體數(shù)值結(jié)果不在此展示。

9 結(jié)論與討論

FWI是勘探地震中地震波反演成像的核心理念。它將地下介質(zhì)視為背景+散射(點)的模式,并試圖對實測炮集中所有波現(xiàn)象進行基于振幅的波形逼近基礎(chǔ)上的反演,其目標(biāo)是估計全或?qū)挷〝?shù)帶的速度場。但勘探地震中的實踐表明,僅有地表觀測的地震數(shù)據(jù),不對勘探地震介質(zhì)的特征和波場特征加以利用,即使僅用透射波進行基于低頻長偏移距數(shù)據(jù)的背景參數(shù)估計,FWI方法也不很穩(wěn)定(視數(shù)據(jù)情況達到收斂或不收斂)。

我們認(rèn)為,一定假設(shè)下的地震波正演方法可以預(yù)測實際地震波場是地震波反演最基本的要求。既然無論如何選擇地震波正演模擬方法都無法完全預(yù)測實測數(shù)據(jù)中的波場,轉(zhuǎn)而追求盡可能好地模擬其中的部分特征波場、進行基于特征波場的反演成像是合理的邏輯。特征表達參數(shù)模型和特征表達數(shù)據(jù)體得到二者之間聯(lián)系更密切的關(guān)系,可以將強非線性的經(jīng)典FWI問題轉(zhuǎn)化為一個較凸的反問題并進行求解?；谏鲜隼砟钐岢龅奶卣鞑ǚ囱莩上?CWI)方法,其中心工作在于特征波場的提取,包括三方面的含義:①特征波型的提取(波型分解);②壓縮感知框架下特征波場的提取(時空局部的方向波場的分解);③反射子波特征的提取(子波分解)。特征波場與反演參數(shù)之間的關(guān)系更緊密,基于特征波場的反問題凸性更好。另外,特征波成像中梯度的預(yù)條件要相對簡單,因為特征波對應(yīng)的梯度干擾大大減弱。

CWI包括四個關(guān)鍵步驟:①用壓縮感知方法提取所要的特征波場;②基于到達時(相位)的透射波波動理論層析反演;③基于到達時(相位)的一次反射波波動理論層析反演;④主要針對方位角度反射系數(shù)的最小二乘疊前深度偏移成像。到目前為止,參數(shù)的低波數(shù)成分估計方法基本確定,就是透射層析和反射波層析,其中的反演方法框架已較為成熟(具體計算中還會有一些問題),但是對應(yīng)的到達時(相位)差測量是一個核心問題。實質(zhì)上,FWI中的周期跳躍就是因為沒有解決好對應(yīng)同相軸的時差(相位)測量問題,否則,周期跳躍問題是不存在的。時差(相位)測量一定會受到振幅的影響,這也是物理上不可避免的,我們能做的只是減少這種影響。我們將時差測量問題作為CWI的核心問題來對待。

關(guān)于參數(shù)的高波數(shù)成分估計,我們認(rèn)為目前還沒有什么公認(rèn)的方法技術(shù)路線。FWI方法假設(shè)地下介質(zhì)是背景+散射(體),而且在一個迭代過程中,從低波數(shù)成分到高波數(shù)成分全部估計出來。實踐證明,FWI方法中高波數(shù)分量估計的收斂性是十分沒有保障的。即便是經(jīng)典的最小二乘RTM(LS_RTM)方法,也是這樣假設(shè)介質(zhì)的分布,在Born近似下進行數(shù)據(jù)逼近的反演,收斂性也沒有什么保障。LS_RTM經(jīng)過這么多年的研究始終不能大規(guī)模處理實際數(shù)據(jù)也說明了這種技術(shù)的局限性。在射線+Born近似下的角度反射系數(shù)最小二乘反演估計可以算是利用了勘探地震地下介質(zhì)特征的高波數(shù)分量估計方法,但是在小尺度散射(繞射)體分布比較多的探區(qū),這樣的成像方法顯然是不合適的。兼顧反射波和散(繞)射波成像的高波數(shù)成分估計方法目前是有必要發(fā)展的,我們正試圖沿這樣的路線發(fā)展出比較實用的高波數(shù)分量估計方法。

[1] SHIN C,CHA Y H.Waveform inversion in the Laplace-Fourier domain[J].Geophysical Journal International,2009,177(3):1067-1079

[2] SHIN C,MIN D.Waveform inversion using a log arithmic wavefield[J].Geophysics,2006，71(3):R31-R42

[4] WU R S,HU C,WANG B.Nonlinear sensitivity operator and inverse thin-slab propagator for tomographic waveform inversion[J].Expanded Abstracts of 84thAnnual Internat SEG Mtg,2014:928-933

[5] SHEN P,SYMES W W.Automatic velocity analysis via shot profile migration[J].Geophysics,2008,73(5):49-59

[6] VAN LEEUWEN T,HERRMANN F J.Mitigating local minima in full-waveform inversion by expanding the search space[J].Geophysical Journal International,2013,195(1):661-667

[7] VAN LEEUWEN T,HERRMANN F J,PETERS B.A new take on fwi-wavefield reconstruction inversion[R].Amsterdam,Netherlands:76thEAGE Conference and Exhibition,2014

[8] 王華忠,馮波,劉少勇,等.疊前地震數(shù)據(jù)特征波場分解、偏移成像與層析反演[J].地球物理學(xué)報,2015,58(6):2024-2034 WANG H Z,FENG B,LIU S Y,et al.Characteristicwavefield decomposition,imaging and inversion with prestack seismic data[J].Chinese Journal of Geophysics,2015,58(6):2024-2034

[9] 王華忠,馮波,王雄文,等.壓縮感知及其在勘探地震中的應(yīng)用[J].石油物探,2016,55(4):467-474 WANG H Z,FENG B,WANG X W.Compressed sensing and its application in seismic exploration[J].Geophysical Prospecting for Petroleum,2016,55(4):467-474

[10] 馮波.數(shù)據(jù)域波動方程層析速度反演方法研究[D].上海:同濟大學(xué),2015 FENG B.Data domain wave equation tomography for velocity inversion[D].Shanghai:Tongji University,2015

[11] ROCKAFELLAR R T.Lagrange multipliers and optimality[J].SIAM review,1993,35(2):183-238

[12] TARANTOLA A.Linearized inversion of seismic reflection data[J].Geophysical Prospecting,1984,32(6):998-1015

[13] BUNKS C,SALECK F M,ZALESKI S,et al.Multiscale seismic waveform inversion[J].Geophysics,1995,60(5):1457-1473

[14] SIRGUE L,PRATT R G.Efficient waveform inversion and imaging:a strategy for selecting temporal frequencies[J].Geophysics,2004,69(1):231-248

[15] VIRIEUX J,OPERTO S.An overview of full-waveform inversion in exploration geophysics[J].Geophysics,2009,74(6):WCC1-WCC26

[16] CARA M,LéVQUE J J.Waveform inversion using secondary observables[J].Geophysical Research Letters,1987,14(10):1046-1049

[17] LUO Y,SCHUSTER G T.Wave-equation traveltime inversion[J].Geophysics,1991,56(5):645-653

[18] GEE L S,JORDAN T H.Generalized seismological data functionals[J].Geophysical Journal International,1992,111(2):363-390

[19] FICHTNER A,KENNETT B L N,IGEL H,et al.Theoretical background for continental-and global-scale full-waveform inversion in the time-frequency domain[J].Geophysical Journal International,2008,175(2):665-685

[20] VAN LEEUWEN T,MULDER W A.A correlation-based misfit criterion for wave-equation traveltime tomography[J].Geophysical Journal International,2010,182(3):1383-1394

[21] MOGHADDAM P P,MULDER W A.The diagonalator:inverse data space full waveform inversion[J].Expanded Abstracts of 82ndAnnual Internat SEG Mtg,2012:1-6

[22] WARNER M,GUASCH L.Adaptive waveform inversion:FWI without cycle skipping[R].Amsterdam,Netherlands:76thEAGE Conference & Exhibition,2014

[23] DE HOOP M V,VAN DER HILST R D,SHEN P.Wave-equation reflection tomography:annihilators and sensitivity kernels[J].Geophysical Journal International,2006,167(3),1332-1352

[24] WEIBULL W,ARNTSEN B,NILSEN E.Initial velocity models for full waveform inversion[J].Expanded Abstracts of 82ndAnnual Internat SEG Mtg,2012:1-4

[25] TARANTOLA A.Inverse problem theory and methods for model parameter estimation[M].Paris,France:SIAM,2005:20-37

[26] 赫孝良,葛照強編著.最優(yōu)化與最優(yōu)控制[M].西安:西安交通大學(xué)出版社,2015:135-155 HE X L,GE Z Q.Optimization and optimal control[M].Xi’an:Xi’an Jiaotong University Press,2015:135-155

[27] BEDNAR J B,SHIN C,PYUN S.Comparison of waveform inversion,part 2:phase approach[J].Geophysical Prospecting,2007,55(4):465-475

[28] SOUBARAS R,GRATACOS B.Velocity model building by semblance maximization of modulated-shot gathers[J].Geophysics,2007,72(5):U67-U73

[29] WARNER M,GUASCH L.Robust adaptive waveform inversion[J].Expanded Abstract of 85thAnnual Internat SEG Mtg,2015:1059-1063

[30] SUKJOON P,SHIN C,BEDNAR J B.Comparison of waveform inversion:part 2:amplitude approach[J].Geophysical Prospecting,2007,55(4):477-485

[31] KREIMER N,STANTON A,SACCHI M D.Tensor completion based on nuclear norm minimization for 5D seismic data reconstruction[J].Geophysics,2013,78(6):V273-V284

[32] OROPEZA V,SACCHI M D.Simultaneous seismic data denoising and reconstruction via multichannel analysis[J].Geophysics,2011,76(3):V25-V32

[33] RECHT B.A simpler approach to matrix completion[J].Journal of Machine Learning Research,2011,12(4):3413-3430

[34] SYMES W W,VERSTEEG R.Velocity model determination using differential semblance optimization[J].Expanded Abstracts of 63rdAnnual Internat SEG Mtg,1993:696-699

[35] SHEN P.Differential semblance velocity analysis via shot profile migration[J].Expanded Abstracts of 75thAnnual Internat SEG Mtg,2005:2249-2253

[36] 劉玉金.復(fù)雜介質(zhì)條件下基于反演的地震成像方法研究[D].青島:中國石油大學(xué)(華東),2014 LIU Y J.Seismic imaging based on inversion in complex media[D].Qingdao:China University of Petroleum (East China),2014

[37] LIU Y,SYMES W W,LI Z.Extended reflection waveform inversion via differential semblance optimization[J].Expanded Abstracts of 84thAnnual Internat SEG Mtg,2014:1149

[38] BAKER P M,GERRITSEN S,CAO Q,et al.3D RTM-based wave-path tomography:theory and applications to synthetic and field data[J].Expanded Abstracts of 85thAnnual Internat SEG Mtg,2015:5259-5264

[39] BAKER P M,GERRITSENS,CAO Q,et al.3D RTM-based wave-path tomography tested at a realistic scale[R].Madrid,Spain:77thEAGE Conference & Exhibition-Workshops,2015

[40] SHAN G,WANG Y.RTM based wave equation migration velocity analysis[J].Expanded Abstracts of 83rdAnnual Internat SEG Mtg,2013:4726-4731

[41] YANG T,SHRAGGE J,SAVA P.Illumination compensation for image-domain wavefield tomography[J].Geophysics,2013,78(5):U65-U76

[42] FEI W,WILLIAMSON P.On the gradient artifacts in migration velocity analysis based on differential semblance optimization[J].Expanded Abstracts of 80thAnnual Internat SEG Mtg,2010:4071-4076

(編輯：戴春秋)

The theoretical framework of characteristic wave inversion imaging

WANG Huazhong,FENG Bo,WANG Xiongwen,HU Jiangtao,LIU Shaoyong,LI Hui,ZHOU Yang

(WavePhenomenaandInversionImagingGroup(WPI),Tongjiuniversity,Shanghai200092,China)

The core of seismic wave inversion imaging is to translate solving a strong nonlinear inverse problem to a further convex inverse problem.The reason of strong nonlinearity in inverse problem lies in the far from linear relation between the real data and the model parameters and the forward problem cannot accurately predict the full wave field.Therefore,we put forward the theoretical framework of characteristic wave inversion (CWI).Its basic idea is not to seek predicting the whole wave phenomena in the real data,but to aspire for interpreting characteristic wave field.Comparing with the full waveform inversion (FWI),CWI is consisted of four basic steps:characteristic wave field (CWF) extracting;wave equation based transmitted wave tomography;wave equation based reflection wave tomography;least square prestack depth migration.CWI is an inversion-imaging framework based on the characteristic wave phenomena in the real data.CWF extraction is the key step,which has three levels meaning:decomposing of the wave phenomena (such as vector wave filed is decomposed into scalar wave field,primaries and multiples are decomposed);decomposing of the temporally and spatially local,single seismic phase,directional and bandlimited wave field from real data;decomposing of the wavelet on the events into the travletime,phase and something else.The compressive sensing method is used in CWF extraction.Wave equation based transmitted wave tomography is based on the first arrival and early arrival traveltime.The inversion imaging based on the disturbance wave field firstly deal with the reflection wave from the layered strata,and then the scattered wave and diffracted wave.The background velocity estimation in the data domain mainly depends on the traveltime and phase,and the effect of the amplitude of the wavelet on the traveltime and phase measuring is removed as possible.The traveltime or phase information can be only used in the background velocity estimation in the imaging domain;it is not suitable for using the amplitude information.For the parameter disturbance estimation,the azimuth/opening-angle reflectivity estimation should be firstly implemented,and then the scattering strength is estimated.Characteristic Wave Inversion Imaging is a realistic procedure to push the classical FWI into the practical application.The numerical experimentation proves the statements.

inversion imaging,characteristic wave field,characteristic wave inversion,nonlinear problem,convex problem,data domain inversion,imaging domain inversion

2016-10-08;改回日期：2016-12-05。

王華忠(1964—),男,教授,主要從事地震波傳播、地震數(shù)據(jù)分析和地震波反演成像方面的研究工作。

國家自然科學(xué)基金項目(41374117)和國家油氣重大專項(2011ZX05005-005-008HZ,2011ZX05006-002,2011ZX05023)聯(lián)合資助。

P631

A

1000-1441(2017)01-0038-12

10.3969/j.issn.1000-1441.2017.01.005

This research is financially supported by the National Natural Science Foundation of China (Grant No.41374117) and the National Science and Technology Major Project of China (Grant Nos.2011ZX05005-005-008 HZ,2011ZX05006-002,2011ZX05023).