吳成梁,王華忠,胡江濤,馬建波

(1.波現(xiàn)象與智能反演成像研究組(WPI),同濟(jì)大學(xué)海洋與地球科學(xué)學(xué)院,上海200092;2.油氣藏與開發(fā)工程國家重點實驗室,成都理工大學(xué),四川成都610059;3.中國石油化工股份有限公司中原油田分公司物探研究院,河南濮陽457000)

隨著“兩寬一高”地震數(shù)據(jù)采集技術(shù)的發(fā)展和進(jìn)步,油氣地震勘探逐步進(jìn)入了高分辨率、高保真的反演成像階段,同時對偏移成像技術(shù)要求也越來越高。Bayes估計理論提供了地震波反演成像的基本思想。無論是去噪、偏移成像,還是層析建模等方法技術(shù),都可在(也應(yīng)該在)Bayes反演框架下討論,因為地震數(shù)據(jù)總受隨機(jī)噪聲的干擾,且實測數(shù)據(jù)總是不完全。

根據(jù)Bayes估計理論,假設(shè)給定的正演算子能夠預(yù)測觀測數(shù)據(jù),在預(yù)測噪聲和反演解具有先驗統(tǒng)計分布特征(一般為高斯分布或廣義高斯分布)基礎(chǔ)上,定義一定范數(shù)意義下的預(yù)測誤差,尋找預(yù)測誤差最小的反演解。根據(jù)統(tǒng)計分布特征,一般數(shù)據(jù)誤差采用L2范數(shù),模型約束采用L1或Cauchy范數(shù)等。若預(yù)測噪聲為高斯分布,此時的Bayes估計可以在最小二乘意義下實現(xiàn);進(jìn)一步地,根據(jù)正算子的線性與否,選擇用線性最小二乘方法或非線性最小二乘方法求解反演問題。

反問題的本質(zhì)是Bayes框架下的后驗概率密度最大化的決策或估計問題。更具實際計算可行性的使代價函數(shù)最小的優(yōu)化問題的解是反問題的本質(zhì)。再進(jìn)一步可落實到算子投影上,TARANTOLA[1]定義了參數(shù)空間及其對偶空間、數(shù)據(jù)空間及其對偶空間,原空間與對偶空間由協(xié)方差算子聯(lián)系在一起,由伴隨算子實現(xiàn)數(shù)據(jù)空間向模型空間的映射,得到模型參數(shù)的估計結(jié)果。從原空間到對偶空間的投影過程為數(shù)據(jù)和模型的特征表達(dá)過程。模型協(xié)方差算子包含地下介質(zhì)的結(jié)構(gòu)特征信息,體現(xiàn)了模型不同點之間的聯(lián)系。數(shù)據(jù)協(xié)方差算子包含了觀測數(shù)據(jù)中相干信號的統(tǒng)計特征,高維數(shù)據(jù)的投影過程等價于數(shù)據(jù)的特征提取過程,通過對數(shù)據(jù)進(jìn)行預(yù)處理,可以提取數(shù)據(jù)中的相關(guān)特征。然而在實際應(yīng)用中很難得到關(guān)于數(shù)據(jù)和模型的協(xié)方差(逆)算子,也沒有具體的解析公式,如何合理利用數(shù)據(jù)和模型的協(xié)方差(逆)算子成為了研究的重點。

在Bayes反演框架下,TARANTOLA[1]全面深入地分析了地震波全波形反演(FWI)的理論問題。理論上,FWI能夠獲得地下全波數(shù)帶參數(shù)場分布。然而經(jīng)典的FWI是強(qiáng)非線性反演問題,實際地震數(shù)據(jù)的地震波FWI受多種因素(比如:缺失低頻長偏移距數(shù)據(jù)、子波未知、正演算子無法模擬所有波現(xiàn)象、數(shù)據(jù)噪聲類型未知等)的制約,因而會出現(xiàn)收斂到局部極值或者不收斂的情況。解決此類問題的基本辦法是將這個強(qiáng)非線性問題轉(zhuǎn)化為一系列線性問題或一個較凸的問題求解[2]。最小二乘逆時偏移(LS_RTM)[3-5]是FWI中典型的線性化子問題,是在背景速度固定的情況下,根據(jù)模擬數(shù)據(jù)與觀測數(shù)據(jù)的波形差異估計地下介質(zhì)的高波數(shù)擾動。從反演的角度理解,常規(guī)的疊前深度偏移成像(PSDM)可以認(rèn)為是最小二乘偏移成像的第1次迭代結(jié)果。在FWI和LS_RTM中,學(xué)者們對模型參數(shù)的正則化[6-7]進(jìn)行了深入探討,提出了很多有效的方法,促進(jìn)了成像質(zhì)量的提高。但是,他們對數(shù)據(jù)正則化的關(guān)注不充分,如何將數(shù)據(jù)質(zhì)量的信息融合在成像過程中以進(jìn)一步提高成像質(zhì)量,是一個值得重視的問題。

本文首先從概率論的角度討論了Bayes框架下的地震波反演成像問題,說明了觀測噪聲為高斯白噪的情況下,Bayes估計可以在最小二乘意義下實現(xiàn);分析了吉洪諾夫(Tikhonov)正則化與Bayes估計的等價性;重點分析了連接原空間與對偶空間的協(xié)方差算子的作用,證明了數(shù)據(jù)特征表達(dá)和模型特征表達(dá)的必要性;從加權(quán)最小二乘誤差泛函角度分析了數(shù)據(jù)協(xié)方差(逆)算子的作用,初步探討了數(shù)據(jù)協(xié)方差(逆)算子在偏移成像中的應(yīng)用;提出了采用傾角掃描和動態(tài)時間規(guī)整算法來估計數(shù)據(jù)加權(quán)系數(shù);最后,理論和實際數(shù)據(jù)的數(shù)值實驗結(jié)果表明數(shù)據(jù)協(xié)方差(逆)算子能夠有效提高成像質(zhì)量。

1 Bayes框架下的地震波反演成像

地震波反演成像建立在Bayes框架下估計理論基礎(chǔ)之上。整個地震數(shù)據(jù)處理與反演成像過程(去噪、偏移成像、速度分析、層析建模等)都可認(rèn)為在此框架下實現(xiàn)[8-9]。

實際地震數(shù)據(jù)是一個隨機(jī)過程,實測到的地震數(shù)據(jù)則認(rèn)為是該隨機(jī)過程的一次具體實現(xiàn)。因此,可以將觀測的地震數(shù)據(jù)視為一組隨機(jī)信號,其具有均值、方差等統(tǒng)計特性,滿足一定的概率分布。同樣地,反演的地下模型參數(shù)也被認(rèn)為是隨機(jī)的。假設(shè)有模型空間M和數(shù)據(jù)空間D,在其定義的聯(lián)合空間M×D上,有同時反映觀測數(shù)據(jù)d以及模型參數(shù)m信息的聯(lián)合先驗概率密度ρ(d,m),同樣有反映聯(lián)系觀測數(shù)據(jù)與模型參數(shù)物理規(guī)律的概率密度Θ(d,m)。利用Bayes公式,可得出聯(lián)合空間M×D上的后驗概率密度σ(d,m):

(1)

其中,k為歸一化常數(shù),μ(d,m)為聯(lián)合空間中關(guān)于觀測數(shù)據(jù)與模型參數(shù)的均勻概率密度函數(shù)。均勻概率密度函數(shù)與坐標(biāo)系統(tǒng)有關(guān),表示特定數(shù)據(jù)和模型參數(shù)化下(選定的坐標(biāo)系統(tǒng))概率空間中單位區(qū)域內(nèi)的概率密度函數(shù),代表了對模型空間與數(shù)據(jù)空間的參數(shù)化認(rèn)識。

進(jìn)一步地,關(guān)于模型空間的后驗概率密度σM(m)以及關(guān)于數(shù)據(jù)空間的后驗概率密度σD(d)可由邊緣概率密度公式計算得到:

假設(shè)數(shù)據(jù)空間與模型空間為兩個獨立空間,聯(lián)合空間中的均勻概率密度μ(d,m)為:

(4)

其中,μD(d)和μM(m)分別表示關(guān)于觀測數(shù)據(jù)和模型參數(shù)的均勻概率密度函數(shù)。(4)式表示模型空間中的模型參數(shù)化與數(shù)據(jù)空間中的數(shù)據(jù)參數(shù)化無關(guān)。

假設(shè)觀測數(shù)據(jù)與模型參數(shù)之間的先驗信息無關(guān),公式(1)中的聯(lián)合先驗概率密度ρ(d,m)變?yōu)?

(5)

其中,ρD(d)和ρM(m)分別表示關(guān)于觀測數(shù)據(jù)和模型參數(shù)的先驗概率密度函數(shù)。

地球物理反演(包括其它基于數(shù)據(jù)的反演)的本質(zhì)都是從數(shù)據(jù)中學(xué)習(xí)到有用的信息并據(jù)此得到客觀認(rèn)識。實測數(shù)據(jù)和要估計參數(shù)之間存在一定的關(guān)系是反演的基礎(chǔ),若二者沒有關(guān)系或關(guān)系極弱,則不可能從實測數(shù)據(jù)中得到要估計的參數(shù)信息。據(jù)此,聯(lián)合概率密度Θ(d,m)可用關(guān)于模型參數(shù)的均勻概率密度函數(shù)及其對應(yīng)的條件概率密度函數(shù)來表示:

(6)

其中,θ(d|m)為在給定的模型參數(shù)m的條件下能夠預(yù)測觀測數(shù)據(jù)d的條件概率密度函數(shù),表示要反演的參數(shù)與實測數(shù)據(jù)之間存在的關(guān)系。均勻概率密度函數(shù)μM(m)的引入表示模型參數(shù)的選擇與地震數(shù)據(jù)的理論預(yù)測無關(guān)。

因此,關(guān)于模型空間的后驗概率密度σM(m)變?yōu)?

(7)

(8)

此時的估計結(jié)果是用后驗概率密度函數(shù)得到的加權(quán)平均結(jié)果。只有當(dāng)后驗概率密度函數(shù)是單峰的,譬如高斯或廣義高斯,此期望值才有意義。

(9)

方差越大,參數(shù)估計的精度在概率意義上越低。

顯然,從概率論的角度看,利用公式(8)進(jìn)行參數(shù)估計是合理的。但是目前缺乏有效的方法計算高維模型參數(shù)的聯(lián)合后驗概率密度函數(shù)。(8)式定義的積分也需要巨大的計算量。為此,將求后驗概率密度函數(shù)的期望值轉(zhuǎn)為求取使得后驗概率密度最大化的模型參數(shù)估計結(jié)果作為反演解:

(10)

由于公式(7)中的后驗概率密度函數(shù)較為復(fù)雜,無法直接求解該式對應(yīng)的最優(yōu)化問題。需要進(jìn)一步引入先驗信息。

首先假設(shè)模型空間M和數(shù)據(jù)空間D均為線性空間,而且模型參數(shù)和觀測數(shù)據(jù)在M和D中的采樣沒有先驗認(rèn)識,則均勻概率密度函數(shù)為:

(11)

假設(shè)聯(lián)系觀測數(shù)據(jù)與模型參數(shù)之間的物理規(guī)律不確定性可以忽略,則:

(12)

假設(shè)關(guān)于模型參數(shù)m的先驗信息是具有高斯分布的概率密度函數(shù),則:

(13)

假設(shè)數(shù)據(jù)空間中數(shù)據(jù)d與實測數(shù)據(jù)dobs之差滿足高斯分布:

(14)

融合上述假設(shè)后,模型空間的后驗概率密度σM(m)變?yōu)?

其中,S(m)為誤差泛函。由單調(diào)性可知,后驗概率密度函數(shù)最大化等價于誤差泛函最小化:

(16)

進(jìn)一步地,誤差泛函有以下形式:

(17)

基于Bayes估計理論可以從概率統(tǒng)計的角度理解反演成像,對反演結(jié)果進(jìn)行精確評價。數(shù)據(jù)協(xié)方差(逆)陣和模型協(xié)方差(逆)陣攜帶了數(shù)據(jù)和模型參數(shù)的不確定性信息,據(jù)此可以更方便地進(jìn)行正則化處理。另外,只有當(dāng)觀測數(shù)據(jù)與模型參數(shù)之間存在弱非線性關(guān)系,后驗概率密度函數(shù)為高斯或廣義高斯分布時,后驗概率密度最大化以及據(jù)此發(fā)展出的代價函數(shù)最小化的方法才有可能用于實際計算。線性的或非線性的最小二乘方法是最基本的算法。

2 數(shù)據(jù)協(xié)方差與模型協(xié)方差

從以上分析可以看出,(15)式和(17)式中的數(shù)據(jù)協(xié)方差CD與模型協(xié)方差CM在模型參數(shù)的優(yōu)化求解中具有特別重要的意義,關(guān)于數(shù)據(jù)和模型參數(shù)的先驗信息很多都蘊含在這兩個協(xié)方差算子中。所謂關(guān)于數(shù)據(jù)的正則化和關(guān)于模型參數(shù)的正則化處理主要在于如何利用好這兩個協(xié)方差矩陣。

首先,從算子投影的角度,引入數(shù)據(jù)和模型參數(shù)的原空間以及它們各自的對偶空間。假設(shè)M和D為線性空間,M*和D*分別為其對應(yīng)的對偶空間,m∈M,d∈D,m*∈M*,d*∈D*。令L為從模型空間映射到數(shù)據(jù)空間的線性算子,即:

(18)

LT定義為從D*映射到M*的線性算子,存在如下的對偶內(nèi)積關(guān)系:

(19)

假設(shè)在線性空間Z中存在一個加權(quán)算子W,該算子是一個線性、對稱、正定的映射算子,它將線性空間Z中的元素z映射到其對偶空間Z*中。線性空間Z中任意兩個向量u,z的標(biāo)量積定義為:

(20)

由加權(quán)算子W建立了原空間與對偶空間的映射,其中,W算子的逆算子C=W-1稱為協(xié)方差算子。令z與z*分別為原空間Z與對偶空間Z*的兩個元素,則有:

(21)

因此,可以說協(xié)方差矩陣將一個線性空間投影到其對偶空間中。而轉(zhuǎn)置算子定義了從數(shù)據(jù)對偶空間向模型參數(shù)對偶空間的映射關(guān)系。標(biāo)量積的定義引入了協(xié)方差算子,建立起原空間的元素到對偶空間元素的映射關(guān)系,反之亦然。

從映射的角度看線性反演,核心是建立起原數(shù)據(jù)空間向原模型參數(shù)空間的映射關(guān)系。這需要引入伴隨算子的概念。L的伴隨算子L*定義為從原數(shù)據(jù)空間D映射到原模型參數(shù)空間M的線性算子,且對所有的元素m∈M和d∈D滿足:

(22)

經(jīng)過標(biāo)量積和轉(zhuǎn)置算子的運算可得到:

(23)

因此伴隨算子為:

(24)

最終從原數(shù)據(jù)空間D映射到原模型參數(shù)空間M的映射過程(圖1)為:

(25)

據(jù)此可以看到,數(shù)據(jù)和模型的協(xié)方差(逆)算子連接起了原空間與對偶空間的橋梁。

接著從線性反演的角度來進(jìn)一步理解協(xié)方差算子。假設(shè)聯(lián)系數(shù)據(jù)和模型之間的正問題是線性、確定性的,此時正問題(12)式變成公式(18),則對應(yīng)的Bayes估計后驗概率密度為:

圖1 原空間與對偶空間關(guān)系示意圖解

(26)

誤差泛函為:

(27)

對誤差泛函求導(dǎo):

(28)

對應(yīng)的法方程為:

(29)

法方程的解為:

(30)

已知協(xié)方差陣是對稱正定陣,可以對其進(jìn)行正交分解。數(shù)據(jù)協(xié)方差陣對數(shù)據(jù)向量進(jìn)行變換可以認(rèn)為是對數(shù)據(jù)進(jìn)行了特征表達(dá)。模型參數(shù)對應(yīng)的協(xié)方差陣具有同樣的作用。因此,可以認(rèn)為對偶空間是特征數(shù)據(jù)和特征模型參數(shù)各自形成的線性空間。特征數(shù)據(jù)和特征模型參數(shù)之間的關(guān)系更緊密、或者更線性、或更弱非線性。轉(zhuǎn)置算子建立起它們之間一對一的投影關(guān)系。這對提高反演結(jié)果的穩(wěn)定性非常有益。這是反演問題中正則化處理最本質(zhì)的內(nèi)涵。

分別對模型和數(shù)據(jù)的協(xié)方差陣做正交分解:

(31)

其中,U和V分別為左奇異向量和右奇異向量,Σ為奇異值組成的對角矩陣。

整理(31)式,得到特征表達(dá)后的解為:

(32)

分析(32)式,得到以下認(rèn)識:

1)Vd算子是將數(shù)據(jù)向量d投影到正交坐標(biāo)系中進(jìn)行特征提取,相當(dāng)于對數(shù)據(jù)進(jìn)行預(yù)處理,提取數(shù)據(jù)中的相干信號;在高維數(shù)據(jù)中,可認(rèn)為是對數(shù)據(jù)進(jìn)行方向譜分解,等價于提取局部平面波過程。

2)VL算子是對傳播算子進(jìn)行相應(yīng)處理,VL對應(yīng)局部平面波傳播算子。(VL)T是對預(yù)條件處理后的局部平面波殘差(或數(shù)據(jù))進(jìn)行反投影的算子。

3)Umprior算子是將模型參數(shù)向量mprior投影到對應(yīng)的正交坐標(biāo)系中,等價于在模型參數(shù)場中提取線性結(jié)構(gòu)信息。

對數(shù)據(jù)實施預(yù)條件處理,提出特征波[2]的概念,同時配合相應(yīng)的模型參數(shù)特征表達(dá),是反演成像中更為根本的正則化思想。據(jù)此可以將一個強(qiáng)非線性問題提成凸問題,或?qū)⒁粋€很病態(tài)問題提成臨近良態(tài)問題,降低反問題求解的難度,提高反演成像的精度。

3 數(shù)據(jù)協(xié)方差在反演(偏移)成像中的應(yīng)用

在地震波反演成像中,人們對模型參數(shù)的正則化進(jìn)行了深入的研究,提出了很多有效的方法,但是對數(shù)據(jù)的正則化考慮較少。根據(jù)前面討論的基本概念,我們引入數(shù)據(jù)正則化方法,提出了自適應(yīng)數(shù)據(jù)變化的加權(quán)疊前深度偏移成像方法。

3.1 加權(quán)最小二乘偏移

忽略模型的正則化項,關(guān)于數(shù)據(jù)匹配項的最小二乘誤差泛函可寫成:

(33)

其中,e表示觀測數(shù)據(jù)與模擬數(shù)據(jù)的誤差(殘差)向量。在Bayes反演框架下,誤差向量(寫成分量形式為:e=(e1,e2,…,eN)為隨機(jī)變量,該隨機(jī)變量的協(xié)方差矩陣可表示為:

(34)

當(dāng)誤差變量ei(i=1,2,…,N)各分量之間高斯統(tǒng)計無關(guān)時,Cee變成Σ,更理想地,變成σ2I。(33)式考慮的就是σ2=1時高斯白噪聲情況下的反演成像。根據(jù)線性高斯情況下Bayes估計理論假設(shè),對于一般的正問題預(yù)測及實際噪聲情況,可認(rèn)為Cee滿足經(jīng)典的正態(tài)分布。當(dāng)反演解趨于(等于)真解時,數(shù)據(jù)殘差滿足高斯白噪分布。此時,誤差向量e中沒有任何有效的信號成分。若正問題不能模擬所有信號成分,或噪聲不符合高斯白噪,導(dǎo)致Cee偏離理想的期望情況時,應(yīng)引入?yún)f(xié)方差逆算子懲罰掉數(shù)據(jù)中的相關(guān)信號成分以滿足上述誤差為高斯白噪的假設(shè)。因此需要在地震波成像中將(33)式修改為加權(quán)最小二乘定義的誤差泛函。

加權(quán)最小二乘誤差泛函可表示為:

(35)

其中,Wd為數(shù)據(jù)加權(quán)函數(shù),表示對數(shù)據(jù)殘差中的不同成分懲罰程度不同,體現(xiàn)了數(shù)據(jù)殘差中不同成分的可靠性。依據(jù)上述討論,最理想的加權(quán)函數(shù)應(yīng)取為誤差向量e的協(xié)方差矩陣的逆矩陣:

(36)

此時,加權(quán)最小二乘偏移成像對應(yīng)的法方程的解為:

(37)

對應(yīng)的迭代公式可寫為:

(38)

其中,k為迭代次數(shù),μk為迭代步長。

然而在實際應(yīng)用中,很難直接通過數(shù)據(jù)計算出完整的協(xié)方差陣及其逆矩陣。數(shù)據(jù)協(xié)方差(逆)陣是與觀測系統(tǒng)有關(guān)的矩陣,一般需要根據(jù)其物理含義賦予不同的假設(shè)形式。

3.2 Beam偏移

(39)

若假設(shè)加權(quán)系數(shù)為對角陣:

(40)

此時數(shù)據(jù)殘差中各個分量之間不相關(guān)且權(quán)重不同,權(quán)系數(shù)依據(jù)數(shù)據(jù)中的噪聲水平或數(shù)據(jù)的可靠性選擇合適的權(quán)重。對于噪聲強(qiáng)的數(shù)據(jù)給予小的權(quán)重,可靠性強(qiáng)的數(shù)據(jù),給予大的權(quán)重。

鑒于勘探地震中殘差向量存在空間相關(guān)性(這是地震波場的空間連續(xù)性決定的),可以用殘差向量相鄰道的數(shù)據(jù)相關(guān)性來構(gòu)造加權(quán)函數(shù)。直接在時間空間域考察(殘差)數(shù)據(jù)的空間相關(guān)性比較困難,一般地,需要在變換域按照特征數(shù)據(jù)來確定合適的權(quán)系數(shù)。

根據(jù)地震波在地下介質(zhì)中的傳播特征,可以選擇用波束合成(Beam Forming)或線性Radon變換的方式構(gòu)建權(quán)函數(shù)。該思想利用局部平面波波前到達(dá)不同檢波器的時差,對觀測數(shù)據(jù)進(jìn)行時移疊加,使局部平面波相干疊加而其它信號相互抵消,從而去除地震數(shù)據(jù)中的某些噪聲。這就是勘探地震中Beam偏移的思想,可看作將數(shù)據(jù)向量d投影到正交坐標(biāo)系中進(jìn)行特征提取的過程。比如Gaussian Beam偏移[10-12]是在中心射線的傍軸近似展開,采用局部高斯函數(shù)計算方式改善焦散、陰影區(qū)和臨界區(qū)振幅計算不準(zhǔn)確問題。SUN等[13]采用在局部地震道內(nèi)進(jìn)行束疊加后再偏移的方法,提高了成像結(jié)果的信噪比和計算效率。

實測地震記錄中若存在波動方程無法描述的波現(xiàn)象,將會導(dǎo)致預(yù)測誤差不符合高斯分布假設(shè)。一種有效的做法是將正演算子所能夠描述的波現(xiàn)象從觀測記錄中剝離出來,然后用分離后的波現(xiàn)象進(jìn)行參數(shù)反演[9]。

3.3 自適應(yīng)數(shù)據(jù)加權(quán)偏移

對比公式(37)與公式(39)可以看出,數(shù)據(jù)協(xié)方差矩陣主要作用在常規(guī)偏移成像結(jié)果和Hessian矩陣中。LS_RTM是在Bayes反演框架下一種有效的估計模型參數(shù)高波數(shù)擾動的保真成像技術(shù)。理想的LS_RTM需要估計精確的Hessian逆算子并作用到常規(guī)的成像結(jié)果上。然而由于Hessian矩陣過于龐大,無法直接估計精確的Hessian逆算子。即使可以采用各種迭代法近似求取Hessian逆算子,仍受限于觀測數(shù)據(jù)不完備、正算子不準(zhǔn)確以及地震子波未知等多種因素影響,需要付出巨大的計算代價才能收斂,甚至出現(xiàn)不收斂現(xiàn)象。因此本文忽略Hessian矩陣,僅在常規(guī)的疊前深度偏移成像中考慮應(yīng)用數(shù)據(jù)協(xié)方差。此時加權(quán)最小二乘法方程的解退化為:

(41)

疊前深度偏移成像中誤差向量e是觀測的地表數(shù)據(jù)減去背景波場(直達(dá)波、透射波等)后的一次波(及散射波)。其中包含相干的波場成分,不是隨機(jī)噪聲,更不是白噪聲。此時引入加權(quán)函數(shù)提高成像質(zhì)量的物理解釋是:We將數(shù)據(jù)進(jìn)行了特征表達(dá),各特征成分之間不相關(guān),高斯假設(shè)下相互獨立,滿足了線性最小二乘反演對殘差向量的假設(shè),可以得到無偏和方差最小的反演解,這是高質(zhì)量反演成像的目標(biāo)。

考慮到數(shù)據(jù)協(xié)方差矩陣逆算子的作用,直接用Cee求逆形成加權(quán)系數(shù)矩陣不可行。如何合理選擇加權(quán)算子W是重要的研究內(nèi)容。利用數(shù)據(jù)中蘊含的空間相關(guān)信息近似地確定加權(quán)函數(shù)W是基于數(shù)據(jù)自適應(yīng)加權(quán)的PSDM方法要處理的核心問題。最基本的思想是測量數(shù)據(jù)中蘊含的相關(guān)信息,據(jù)此構(gòu)造加權(quán)算子W。

基于上述理論分析,假設(shè)加權(quán)系數(shù)矩陣為對角矩陣,加權(quán)系數(shù)可由觀測波場中的相干成分近似計算。采用傾角掃描和動態(tài)時間規(guī)整算法確定加權(quán)系數(shù),其基本流程為:

1) 在局部的共炮(CS)/共中心點(CMP)道集中提取局部波場;

2) 在選擇的局部波場中,采用傾角掃描的方式獲取不同道之間的線性相位擾動;

3) 采用動態(tài)時間規(guī)整算法校正非線性相位擾動;

4) 在拉平后的道集上利用相似系數(shù)譜方法估計信息中的相干成分,從而估計出加權(quán)系數(shù);

5) 重復(fù)步驟1)到步驟4),直到估計出完整的加權(quán)系數(shù);

6) 采用公式(41)進(jìn)行偏移成像。

關(guān)于具體的傾角掃描和動態(tài)時間規(guī)整算法可參考文獻(xiàn)[14]。

加權(quán)算子的施加方式有:①直接改造疊前數(shù)據(jù);②作用于Kirchhoff積分公式中的加權(quán)系數(shù)上;③改造波場外推成像方法中的相關(guān)成像條件。本文僅實現(xiàn)了第一種方法,后兩種方法在今后的研究中逐步展開。

3.4 數(shù)值實驗

3.4.1 模擬數(shù)據(jù)

采用圖2所示單層理論模型模擬數(shù)據(jù)對本文方法進(jìn)行測試。圖2a為采用有限差分法模擬的單炮記錄,共601道,道間隔10m,時間采樣點數(shù)為5000,時間采樣率為0.5ms;圖2b為加入隨機(jī)噪聲后的單炮記錄,其中信噪比為2;圖2c為采用本文方法計算的加權(quán)系數(shù),由圖2c可以看到,有效信號(相干信號)存在的地方,加權(quán)系數(shù)接近于1,而在非相干信號區(qū)域,加權(quán)系數(shù)基本為0;圖2d為采用數(shù)據(jù)加權(quán)處理后得到的單炮記錄,由圖2d可以看出,噪聲得到了有效消除。圖3a和圖3b分別給出了數(shù)據(jù)加權(quán)處理前、后的含噪單炮記錄局部放大結(jié)果。由圖3可以看出,原始單炮記錄中包含能量較強(qiáng)的隨機(jī)干擾噪聲,經(jīng)過數(shù)據(jù)加權(quán)處理后,隨機(jī)干擾噪聲得到了消除,較為完整地保留了有效信號。

圖2 單層理論模型模擬數(shù)據(jù)測試結(jié)果a 原始單炮記錄; b 含噪單炮記錄; c 采用本文方法計算的加權(quán)系數(shù); d 數(shù)據(jù)加權(quán)處理后的單炮記錄

3.4.2 實際地震數(shù)據(jù)

采用某工區(qū)二維炮集數(shù)據(jù)測試本文方法的有效性。圖4a為原始單炮記錄,可以看到,原始單炮記錄中含有較多的相干噪聲,深層同相軸被噪聲掩蓋,無法識別;采用本文方法處理后的單炮記錄如圖4b所示,可以看出,噪聲得到有效去除,無論是淺層還是深層,同相軸連續(xù)性更好。對應(yīng)的偏移成像結(jié)果如圖5所示。從圖5可以清晰地看到,采用本文方法處理后的偏移剖面成像質(zhì)量得到明顯提升,同相軸更加清晰連續(xù),成像剖面上噪聲明顯減弱。

采用某工區(qū)三維CMP道集數(shù)據(jù)測試本文方法的效果。圖6a為原始的CMP道集;圖6b為數(shù)據(jù)加權(quán)處理后的CMP道集。由圖6可以看出,淺層能量強(qiáng)的同相軸表現(xiàn)出較強(qiáng)的相似性,在處理后的道集中得到增強(qiáng),而在沒有相干信號的地方,原始道集中含有較多的噪聲,處理后的道集中噪聲得到有效去除。圖7為某條測線采用本文方法處理前、后的偏移成像結(jié)果。對比圖7a和圖7b發(fā)現(xiàn),采用本文方法處理后,成像質(zhì)量得到明顯改善。

圖4 采用本文方法處理前(a)、后(b)的單炮記錄

圖5 采用本文方法處理前(a)、后(b)的偏移成像剖面

圖6 數(shù)據(jù)加權(quán)處理前(a)、后(b)的CMP道集

圖7 某測線采用本文方法處理前(a)、后(b)的偏移成像結(jié)果

4 結(jié)論與討論

“兩寬一高”地震勘探的目的是進(jìn)行更精細(xì)的油藏描述?！皟蓪捯桓摺钡卣饠?shù)據(jù)采集提供了更完善的地震數(shù)據(jù),但是地震波成像技術(shù)的發(fā)展相對滯后,不能充分挖掘出所采集的巨量數(shù)據(jù)中蘊含的與儲層有關(guān)的信息。地震波成像應(yīng)該逐漸地奠定在Bayes估計理論基礎(chǔ)上,從信息綜合的角度提高地震波成像的質(zhì)量。FWI和LS_RTM充分關(guān)注了模型參數(shù)的正則化,促進(jìn)了成像質(zhì)量的提高。但是,對數(shù)據(jù)正則化的關(guān)注有所欠缺,將數(shù)據(jù)中蘊含的先驗信息融合在成像過程中進(jìn)一步提高成像質(zhì)量是值得重視的。

為此,本文首先在Bayes框架下討論了地震波反演成像問題。在地震波正問題(波場和數(shù)據(jù)間的關(guān)系)是弱線性的假設(shè)下,反演成像的核心問題是基于后驗概率密度最大化原則的參數(shù)估計問題。若觀測噪聲(嚴(yán)格地講是預(yù)測誤差)為高斯白噪,模型參數(shù)的概率分布為高斯分布情況下,最小二乘估計可以在最小二乘意義下實現(xiàn)。但是,實際數(shù)據(jù)并不滿足上述假設(shè),即預(yù)測誤差滿足有色高斯分布、甚至不滿足高斯分布。此時,引入數(shù)據(jù)協(xié)方差(逆)矩陣是提高反演成像質(zhì)量所必須的。數(shù)據(jù)和模型的協(xié)方差(逆)矩陣連接了各自的原空間與對偶空間,通過協(xié)方差矩陣進(jìn)行數(shù)據(jù)和模型參數(shù)的正交投影,相當(dāng)于分別對數(shù)據(jù)和模型參數(shù)進(jìn)行了特征表達(dá)。由此可以在數(shù)據(jù)預(yù)條件的基礎(chǔ)上提取特征波場,并對模型參數(shù)進(jìn)行特征表達(dá),建立起更為緊密的特征波場與特征參數(shù)之間的聯(lián)系,構(gòu)建更凸的反問題。這是進(jìn)行數(shù)據(jù)自適應(yīng)加權(quán)疊前深度偏移的數(shù)學(xué)物理理論基礎(chǔ)。

本文從最小二乘誤差泛函角度分析了數(shù)據(jù)協(xié)方差(逆)算子的作用,指出在地震波反演成像中加權(quán)最小二乘引入的必要性。在假設(shè)數(shù)據(jù)協(xié)方差矩陣為對角陣的基礎(chǔ)上,針對疊前深度偏移成像,提出了一種自適應(yīng)數(shù)據(jù)加權(quán)的數(shù)據(jù)預(yù)條件處理方法。該方法通過傾角掃描和動態(tài)時間規(guī)整方法自動計算數(shù)據(jù)加權(quán)系數(shù),將其融合在偏移成像過程中。理論和實際資料測試結(jié)果驗證了該方法的有效性。

“兩寬一高”地震勘探是一項綜合性技術(shù),采集技術(shù)還有巨大的發(fā)揮空間,反演成像技術(shù)的差距更大?；贐ayes理論,充分地融合關(guān)于數(shù)據(jù)的先驗信息和關(guān)于模型參數(shù)的先驗信息是提高成像質(zhì)量的必由之路。地震波反演成像中對數(shù)據(jù)不確定引起的成像結(jié)果的不確定性研究還處于初始階段,本文的研究工作僅僅是較為初步的探索。

致謝:感謝中石油勘探開發(fā)研究院及西北分院、中海油研究院和湛江分公司、中國石油化工股份有限公司石油物探技術(shù)研究院和勝利油田分公司對波現(xiàn)象與智能反演成像研究組(WPI)研究工作的資助與支持。