地震波反演成像方法與技術(shù)核心問題分析

王華忠,馮 波,王雄文,胡江濤,李 輝

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

地震波反演成像方法與技術(shù)核心問題分析

王華忠,馮 波,王雄文,胡江濤,李 輝

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

常規(guī)的地震反演成像分為偏移速度分析與層析成像、疊前深度偏移(角度道集產(chǎn)生)和AVA分析/反演3個(gè)重要環(huán)節(jié)，其中的關(guān)鍵技術(shù)是當(dāng)前勘探地震學(xué)中的核心技術(shù)。全波形反演(FWI)是理論意義下十分完善的地震波反演成像理論框架。原則上，F(xiàn)WI可以把上述3項(xiàng)常規(guī)的反演方法技術(shù)合為一體，給出比較理想的反演成像結(jié)果。但是，由于疊前地震數(shù)據(jù)的不完備、地震波正演模擬方法不能很好地模擬實(shí)測(cè)地震波場(chǎng)、初始模型不夠精確、地震子波的未知和空變，使得嚴(yán)格意義下的FWI方法尚不能很好地解決實(shí)際問題。通過疊前地震數(shù)據(jù)和地下介質(zhì)模型的特征表達(dá)，提出把經(jīng)典的FWI分成透射波層析成像、最小二乘疊前深度偏移成像和反射波層析成像3個(gè)線性化反演方法的串聯(lián)，構(gòu)成FWI反演成像的實(shí)用化流程。針對(duì)我國陸上地震數(shù)據(jù)的特點(diǎn)，指出做好淺層速度模型建立、背景速度模型建立、成像道集層析速度模型建立、最小二乘疊前深度偏移成像、張角及界面傾角道集的產(chǎn)生以及小角度成像道集波阻抗反演，是當(dāng)前推進(jìn)FWI反演成像方法技術(shù)應(yīng)用與發(fā)展的關(guān)鍵所在。

反演成像理論;多級(jí)Born近似;全波形反演;逐級(jí)線性化反演;實(shí)用化地震反演成像技術(shù)路線

現(xiàn)代地震勘探中，巖性油氣藏的精細(xì)描述主要是利用高精度采樣(寬方位、寬頻帶和高密度采樣)的高維地震數(shù)據(jù)反演估計(jì)精細(xì)的背景參數(shù)變化(主要是速度、各向異性參數(shù)和Q值)和高保真及高分辨率的高波數(shù)參數(shù)變化(主要是反射系數(shù))，結(jié)合巖石物理知識(shí)，進(jìn)行油氣預(yù)測(cè)。精細(xì)背景參數(shù)變化的估計(jì)主要是利用層析成像技術(shù);高保真及高分辨率的高波數(shù)參數(shù)變化的估計(jì)方法主要是保真疊前深度偏移或最小二乘疊前深度偏移技術(shù)?，F(xiàn)代地震勘探的一個(gè)主要特征是整個(gè)成像處理流程逐步進(jìn)入以地震波反演成像為主的階段，其中的每一個(gè)環(huán)節(jié)基本上都可以抽象為建立正演預(yù)測(cè)模型和Bayes參數(shù)估計(jì)兩個(gè)關(guān)鍵步驟;另一個(gè)主要特征是現(xiàn)代的地震反演成像技術(shù)以巖性油氣藏的描述為目標(biāo)。

一般意義下的地震波成像所指的就是地震波偏移成像。它是在假設(shè)偏移速度場(chǎng)已知的情況下，利用在該速度場(chǎng)中計(jì)算的波場(chǎng)反向外推算子，把地表記錄的波場(chǎng)外推到繞射點(diǎn)上，用適當(dāng)?shù)某上駰l件提取成像值，目標(biāo)是定量地定位反射系數(shù)出現(xiàn)的空間位置，并定性地保持反射系數(shù)的相對(duì)強(qiáng)度，尤其是角度反射系數(shù)的相對(duì)強(qiáng)度。

在相對(duì)保真的疊前成像道集上，進(jìn)行地震子波屬性的提取及基于地震子波屬性的儲(chǔ)層識(shí)別與描述;也可以進(jìn)行AVA分析/反演，進(jìn)而開展半定量的儲(chǔ)層識(shí)別與描述。顯然，有比較合適的疊前地震數(shù)據(jù)、正確的背景速度場(chǎng)和保真的角度成像道集是進(jìn)行后續(xù)的儲(chǔ)層解釋的基礎(chǔ)。

地震波偏移成像與地震波反演成像之間的區(qū)別是明顯的。原則上，前者定量地估計(jì)反射點(diǎn)的位置，定性地估計(jì)角度反射系數(shù)的大小;后者試圖定量地估計(jì)全波數(shù)帶的速度場(chǎng)(多參數(shù)反演情況下也包括其它參數(shù)的全波數(shù)帶成分的估計(jì))。地震波反演成像的代表性技術(shù)是全波形反演(full waveform inversion，F(xiàn)WI)，它的目標(biāo)是估計(jì)全波數(shù)帶的速度場(chǎng)(也可以估計(jì)全波數(shù)帶的其它參數(shù)場(chǎng))。

地震波反演成像的基本思想是Bayes估計(jì)理論。Tarantola詳細(xì)而深入地分析了地震波全波形反演的理論問題[1-2]。Bayes估計(jì)理論更是現(xiàn)代信號(hào)分析中的核心內(nèi)容[3-4]。在整個(gè)地震數(shù)據(jù)處理與反演成像過程中，去噪、反褶積、高維數(shù)據(jù)規(guī)則化、一維波阻抗反演、AVA(疊前)彈性參數(shù)反演、最小二乘疊前偏移成像、層析成像等都是在此框架下進(jìn)行的。根據(jù)Bayes估計(jì)理論，首先假設(shè)存在一個(gè)正演預(yù)測(cè)算子，它可以預(yù)測(cè)觀測(cè)數(shù)據(jù)。若預(yù)測(cè)噪聲為高斯白噪，此時(shí)的Bayes估計(jì)可以在最小二乘意義下實(shí)現(xiàn)，此為最大似然估計(jì)方法。假設(shè)已知關(guān)于模型的先驗(yàn)概率分布，可以建立一個(gè)聯(lián)合的后驗(yàn)概率分布，數(shù)據(jù)誤差一般用L2范數(shù)，模型約束一般用L1或Cauchy范數(shù)等;根據(jù)正算子的線性與否，選擇用線性最小二乘方法或非線性最小二乘方法求解反演問題。此為經(jīng)典的Bayes估計(jì)方法。

預(yù)測(cè)觀測(cè)數(shù)據(jù)的正算子十分重要。褶積算子用于反褶積;一維波動(dòng)方程用于預(yù)測(cè)自激自收的地震道(波阻抗反演);自回歸(auto regressive,AR)模型算子常用于線性信號(hào)預(yù)測(cè)(壓制噪聲);旅行時(shí)計(jì)算算子用于預(yù)測(cè)實(shí)測(cè)數(shù)據(jù)中的旅行時(shí)(旅行時(shí)層析成像);Zoeppritz方程及其各種簡(jiǎn)化形式預(yù)測(cè)角度反射系數(shù)(AVA彈性參數(shù)反演);Fourier變換算子用于預(yù)測(cè)規(guī)則無假頻數(shù)據(jù)的譜(高維數(shù)據(jù)規(guī)則化);動(dòng)校正時(shí)距關(guān)系預(yù)測(cè)實(shí)測(cè)CMP道集中的雙曲時(shí)距關(guān)系(估計(jì)均方根速度);各種形式的波動(dòng)方程預(yù)測(cè)實(shí)際觀測(cè)的炮記錄(全波形反演)?？梢钥闯?，每一類正算子對(duì)應(yīng)一類反演問題。正算子建立了要估計(jì)的參數(shù)與觀測(cè)數(shù)據(jù)之間的線性或非線性關(guān)系。實(shí)際上，很多情況下正演結(jié)果并不能很好地預(yù)測(cè)實(shí)測(cè)數(shù)據(jù)，從而導(dǎo)致反演結(jié)果精度降低，甚至反演過程根本不收斂。

反演問題是典型的不適定問題。在地球物理反演中，解的存在性是有物理保證的，問題主要是解的不唯一性和計(jì)算過程中解的不穩(wěn)定性。二者具有關(guān)聯(lián)性，計(jì)算的不穩(wěn)定性加劇了解的不唯一性，但最根本的還是解的不唯一性問題。反演問題解的不唯一性本質(zhì)上是由數(shù)據(jù)空間向模型空間映射的非線性性引起的，甚至是由模型參數(shù)變化的復(fù)雜性和描述波場(chǎng)與參數(shù)場(chǎng)關(guān)系的正算子決定的。不唯一性只可以減弱，不可能完全消除。反演解的不唯一性和反演解的精度可以認(rèn)為是一個(gè)問題的兩個(gè)方面。減弱反演解的不唯一性的關(guān)鍵是減弱由數(shù)據(jù)空間向模型空間映射的非線性性、提高初始模型的精度、增加關(guān)于模型的先驗(yàn)信息，這些措施可以有效地提高反演解的精度。目前，F(xiàn)WI技術(shù)發(fā)展過程中提出的所有方法都是基于上述3點(diǎn)的。

20世紀(jì)80年代，地震波反演處于熱潮期，很多學(xué)者對(duì)全波形反演理論做了十分扎實(shí)的基礎(chǔ)研究工作[5-10]。但受當(dāng)時(shí)野外地震數(shù)據(jù)觀測(cè)技術(shù)和高性能計(jì)算機(jī)技術(shù)的局限，無法將提出的理論方法轉(zhuǎn)為實(shí)用。此后，一維波阻抗反演、AVA疊前彈性參數(shù)反演這些計(jì)算需求小的反演方法得到了不斷發(fā)展。但從理論上看，這些反演方法都不具備高維疊前反演的潛力，對(duì)于介質(zhì)分布情況的假設(shè)過多。從實(shí)用性上看，它們的反演精度從來沒有被認(rèn)真地考究過。可以說，這些反演結(jié)果并沒有對(duì)地震解釋產(chǎn)生過實(shí)質(zhì)性的貢獻(xiàn)。

近幾年，高維疊前地震數(shù)據(jù)的全波形反演理論與方法重新成為了勘探地球物理領(lǐng)域的核心議題，這首先得益于地震數(shù)據(jù)采集技術(shù)和大規(guī)模高性能計(jì)算機(jī)技術(shù)的進(jìn)步，其次是很多研究者對(duì)地震波高維疊前反演研究的興起做了大量的鋪墊工作[11-16]。寬方位、高密度和寬頻帶(尤其是低頻)數(shù)據(jù)的采集技術(shù)和GPU/CPU異構(gòu)計(jì)算機(jī)群技術(shù)為FWI技術(shù)的可能應(yīng)用奠定了很好的基礎(chǔ)。

但是，經(jīng)典FWI技術(shù)的實(shí)際應(yīng)用遇到了巨大的困難，真正成功的應(yīng)用實(shí)例尚不多見。主要體現(xiàn)在疊前數(shù)據(jù)的不完備、地震波正演模擬方法的不完善、初始模型不夠精確、地震子波的未知和空變等4個(gè)方面。近幾年研究者們逐步認(rèn)識(shí)到經(jīng)典的FWI很難得到合適的反演解，于是提出了各種各樣的FWI變種方法，基本的思想是盡量克服上述4個(gè)方面因素的影響，使得反演結(jié)果更加穩(wěn)定，反演過程的收斂性更強(qiáng)。

從理論發(fā)展的角度看，目前以地震波偏移成像為主的疊前地震數(shù)據(jù)成像處理發(fā)展到以FWI技術(shù)為代表的地震波反演成像是無法回避的。但是，F(xiàn)WI方法技術(shù)的實(shí)現(xiàn)還是要從系統(tǒng)的觀點(diǎn)出發(fā)，依據(jù)實(shí)際數(shù)據(jù)情況構(gòu)建出合適的技術(shù)流程，達(dá)到FWI的目的，即得到全(寬)波數(shù)帶的速度估計(jì)結(jié)果(也可以包括其它參數(shù)的估計(jì)結(jié)果)。理論上講，F(xiàn)WI技術(shù)的實(shí)現(xiàn)等于把疊前數(shù)據(jù)采集、反演成像和儲(chǔ)層解釋有機(jī)地融合在一起。

當(dāng)前對(duì)FWI方法技術(shù)的定位，一般地認(rèn)為它可以提供滿足逆時(shí)深度偏移(RTM)成像要求的更為精細(xì)的速度場(chǎng)。實(shí)質(zhì)上，F(xiàn)WI的真實(shí)目的是得到全(寬)波數(shù)帶的參數(shù)估計(jì)，類似RTM的偏移成像方法是輔助FWI得到波數(shù)帶更寬的參數(shù)估計(jì)結(jié)果，而不是FWI僅僅提供一個(gè)更好的偏移速度給RTM[17]。據(jù)此可以看出，F(xiàn)WI是期望從疊前數(shù)據(jù)中定量地估計(jì)彈性參數(shù)(甚至巖性參數(shù)和流體參數(shù))的技術(shù)，但是能否達(dá)到這樣宏大的目標(biāo)，疊前數(shù)據(jù)的質(zhì)量是關(guān)鍵因素。

因此，本文首先從疊前地震數(shù)據(jù)的抽象認(rèn)識(shí)入手，提出特征波的概念以滿足不同層次的FWI對(duì)數(shù)據(jù)的要求，接著提出對(duì)參數(shù)空間及地震波傳播的抽象認(rèn)識(shí)，建立反演參數(shù)和特征波場(chǎng)之間的關(guān)系，并據(jù)此盡量強(qiáng)化二者之間的線性關(guān)系，同時(shí)說明正問題線性化的必要性及給出比較精確的初始模型的必要性。在此基礎(chǔ)上，給出地震波反演成像的合理的技術(shù)流程，以達(dá)到實(shí)現(xiàn)全(寬)波數(shù)帶的參數(shù)估計(jì)的目的。

我們認(rèn)為，地震波反演成像的目的是通過彈性參數(shù)的估計(jì)來刻畫和描述儲(chǔ)層。這是下一代地震成像技術(shù)的顯著特征，標(biāo)志著勘探地震學(xué)全面進(jìn)入反演階段。

1 油氣地震勘探中的關(guān)鍵問題分析

油氣地震勘探的終極目標(biāo)是盡可能精確地描述含油氣儲(chǔ)層。宏觀地講，這是一個(gè)更大的反問題。所依據(jù)的數(shù)據(jù)包括地震數(shù)據(jù)、測(cè)井?dāng)?shù)據(jù)、地質(zhì)數(shù)據(jù)(包括基本地質(zhì)認(rèn)識(shí))、巖石性質(zhì)數(shù)據(jù)(包括基本巖石性質(zhì)測(cè)定)，反演過程歸結(jié)為一個(gè)綜合判斷和決策過程。當(dāng)前，世界上大的地球物理服務(wù)公司都試圖推出這樣的綜合性的軟件系統(tǒng)，可以簡(jiǎn)稱為采集、處理和解釋一體化系統(tǒng)。但這樣的系統(tǒng)在決策階段存在軟肋，解釋方面也是以構(gòu)造解釋為主，含油氣性的判斷甚至還遠(yuǎn)沒有達(dá)到半定量的程度。在反演理論大框架下，向智能化的地震數(shù)據(jù)采集、處理、解釋和決策一體化方向發(fā)展是沒有異議的發(fā)展方向。目前存在如下一些共識(shí)：①地震數(shù)據(jù)采集向?qū)挿轿弧⒏呙芏群蛯掝l帶數(shù)據(jù)采集發(fā)展;②反演成像向全(寬)波數(shù)帶參數(shù)估計(jì)方向發(fā)展的趨勢(shì)已十分明顯;③地震地質(zhì)解釋的全三維化與智能化正逐步展開。這3點(diǎn)的目的都是指向更精確、更智能化地描述儲(chǔ)層。

基于地震資料描述儲(chǔ)層時(shí)，地震波場(chǎng)中攜帶的哪些信息可以對(duì)儲(chǔ)層的存在和儲(chǔ)層的性質(zhì)有所指示，是油氣地震勘探工作者必須思考和弄清楚的問題。這也是反演理論中的基本概念，即數(shù)據(jù)中一定要包含有要反演參數(shù)的變化信息，而且二者的關(guān)聯(lián)性還要足夠強(qiáng)。

目前，可以確認(rèn)的是地震數(shù)據(jù)同相軸的到達(dá)時(shí)與地震波速度、各向異性參數(shù)是有密切關(guān)系的;地震波波阻抗的變化與反射地震波振幅有密切關(guān)系;吸收衰減參數(shù)(Q值)與反射地震波的振幅有明顯的關(guān)系;小尺度地下介質(zhì)體的散射與地震波振幅、頻率和相位存在明顯的關(guān)系。充分利用上述關(guān)系來解決儲(chǔ)層的描述和刻畫是問題的關(guān)鍵所在。目前的地震波反演成像還遠(yuǎn)沒有把上述關(guān)系利用得很徹底。

疊前地震數(shù)據(jù)是波在地下介質(zhì)中傳播后的累積效應(yīng)，因此，利用疊前地震數(shù)據(jù)反演得到的往往也是介質(zhì)性質(zhì)參數(shù)變化的光滑結(jié)果。對(duì)速度、各向異性參數(shù)和Q值的估計(jì)得到的都是這樣的結(jié)果。從局部平面波傳播的角度看，疊前地震數(shù)據(jù)與地下介質(zhì)參數(shù)之間由局部平面波傳播聯(lián)系起來，方位張角成像道集或方位傾角成像道集是油氣儲(chǔ)層描述最核心的基礎(chǔ)數(shù)據(jù)。可以說，它們是成像域的高維疊前數(shù)據(jù)體，直接攜帶了與儲(chǔ)層及周邊巖性變化有關(guān)的信息。到目前為止，在地震反演成像過程中，產(chǎn)生盡可能保真的方位張角成像道集或方位傾角成像道集是進(jìn)一步描述裂縫儲(chǔ)層、孔洞儲(chǔ)層和其它類型巖性儲(chǔ)層的基礎(chǔ)。對(duì)波阻抗變化、裂隙方向和密度、散射體的估計(jì)都需要在這樣的道集上進(jìn)行，利用成像道集上的振幅變化考察的是參數(shù)的高波數(shù)變化成分。

2 疊前地震數(shù)據(jù)的抽象認(rèn)識(shí)

疊前地震反演成像依賴的是單炮道集表示的地震數(shù)據(jù)，一個(gè)探區(qū)的所有單炮道集可以表示為：

(1)

式中：Xs=(xs,ys)和Xr=(xr,yr)分別代表炮點(diǎn)和檢波點(diǎn)坐標(biāo);Ns代表總炮數(shù);R2代表炮點(diǎn)和檢波點(diǎn)分布區(qū)域。

期望這樣的數(shù)據(jù)是寬方位、高密度和寬頻帶(尤其包含低頻成分)的，最好是炮點(diǎn)與檢波點(diǎn)采樣間隔有大致相同的空間采樣率。原因可以由下列定理清楚地表明。

定理[18]：常速情況下，在一階Born近似成立的條件下，在散射體的遠(yuǎn)場(chǎng)范圍內(nèi)，散射勢(shì)F(r)與它的波數(shù)譜f(s,s0)之間存在如下Fourier變換關(guān)系：

(2)

式中：K=(ω/v)·(s-s0)，s0代表入射波方向，s代表出射波方向；r代表散射體的空間分布范圍。(2)式定義了散射勢(shì)F(r)和散射勢(shì)的譜f(s,s0)之間的Fourier變換關(guān)系。

由(2)式可見，如果得到寬波數(shù)帶的散射勢(shì)譜就可以通過Fourier反變換得到高分辨率的散射勢(shì)的估計(jì)。散射勢(shì)與速度擾動(dòng)密切相關(guān)，很多文獻(xiàn)上都有定義。散射勢(shì)的波數(shù)分布范圍由K=(ω/v)·(s-s0)給出。其中，矢量s-s0定義了散射(或反射)張角，觀測(cè)系統(tǒng)決定了張角的范圍，張角范圍越大，波數(shù)譜越寬，對(duì)應(yīng)地，要求地震數(shù)據(jù)的觀測(cè)角度越大，即要有充分長(zhǎng)的偏移距和充分寬的方位角。其中頻率范圍也要求足夠?qū)?，尤其是低頻成分對(duì)散射勢(shì)的低波數(shù)成分貢獻(xiàn)很大。至于疊前地震數(shù)據(jù)中的高密度采樣，主要是從反假頻的需要而提出的，尤其是對(duì)某些特殊的低視速度噪聲需要特別密的空間采樣。當(dāng)然，增加覆蓋次數(shù)、提高信噪比也是加密采樣的基本考慮。

更進(jìn)一步地，從理論上來說，地震波反演成像要求疊前地震數(shù)據(jù)采集系統(tǒng)對(duì)地下任何一個(gè)繞射點(diǎn)(反射點(diǎn))都有廣角度的、角度間隔均勻的、不產(chǎn)生采樣假頻的照明。同時(shí)期望每個(gè)角度的數(shù)據(jù)中僅僅有高斯白噪聲;還期望各角度之間的子波特征保持一致。

目前，比較流行的經(jīng)典意義下的FWI反演是利用變密度聲波方程，甚至常密度聲波方程進(jìn)行地震波正演。當(dāng)然，也有用彈性波方程進(jìn)行地震波正演來預(yù)測(cè)實(shí)測(cè)炮道集的。無論采用何種波動(dòng)方程都不可能很精確地預(yù)測(cè)出實(shí)際數(shù)據(jù)，其中的影響因素實(shí)在太多。不能預(yù)測(cè)的部分變成數(shù)據(jù)殘差，該殘差不符合Gauss分布，使得問題的非線性性變強(qiáng)，解的非唯一性變強(qiáng)。

問題的解決辦法是強(qiáng)調(diào)對(duì)數(shù)據(jù)的預(yù)處理，而不強(qiáng)調(diào)使用復(fù)雜的正演模擬方程。為此，我們提出了特征波的概念。所謂特征波是指單一震相的時(shí)空局部的波，可以是反射波、繞射波或透射波，也可以是某種波的波形、包絡(luò)、相位、甚至是到達(dá)時(shí)。我們認(rèn)為當(dāng)前的反演成像還是要先針對(duì)單一震相的波現(xiàn)象進(jìn)行。

另一方面，我們把單炮地震數(shù)據(jù)抽象地描述為：將地表(也可以是井中或海底(OBC)等處)觀測(cè)到的地震數(shù)據(jù)視為以旅行時(shí)、炮/檢點(diǎn)坐標(biāo)(或中點(diǎn)半偏移距)為變量的函數(shù)，它的基本特征應(yīng)該是由時(shí)距關(guān)系規(guī)定的、由帶限反射子波表現(xiàn)出的同相軸加上一定的隨機(jī)噪聲。速度和各向異性參數(shù)場(chǎng)的背景部分決定同相軸出現(xiàn)的位置，Q值和波阻抗變化決定反射子波的相對(duì)振幅。觀測(cè)數(shù)據(jù)中子波振幅的絕對(duì)值影響因素甚為復(fù)雜，海上地震勘探還比較容易分析清楚，陸上疊前地震數(shù)據(jù)中絕對(duì)振幅的具體含義及影響因素還有待深入分析。這也說明了基于波形振幅差的一類FWI反演技術(shù)難以得到合理的反演解的基本原因。所有基于振幅的反演方法都存在需要講清楚實(shí)測(cè)數(shù)據(jù)的振幅是什么含義的問題。據(jù)此，我們更進(jìn)一步地提出，首先利用同相軸的到達(dá)時(shí)信息反演背景速度、各向異性參數(shù),也可以引入吸收衰減旅行時(shí)反演背景Q值參數(shù)。在此基礎(chǔ)上，進(jìn)行帶吸收衰減的各向異性疊前深度偏移成像，產(chǎn)生保真的方位張角成像道集或方位傾角成像道集，再進(jìn)一步描述和刻畫儲(chǔ)層。這是目前比較現(xiàn)實(shí)的地震波反演成像技術(shù)路線。

更進(jìn)一步地，把單炮道集中的波現(xiàn)象分為透射波、反射波和繞射波等。當(dāng)然也包括多次波和AVA現(xiàn)象。對(duì)應(yīng)于近偏移距的折射波和直達(dá)波的這部分透射波，主要用來進(jìn)行淺層近地表速度層析反演;對(duì)應(yīng)于中長(zhǎng)偏移距的來自中深層的折射波和潛波(Diving Wave)的這部分透射波，主要用來進(jìn)行中深層背景速度的估計(jì)。因?yàn)橥干洳ㄊ遣皇苄〕叨人俣犬惓ｓw的影響的，利用它也只能反演大尺度的背景速度場(chǎng)。目前疊前偏移成像主要用的是一次反射波和/或繞射波。RTM原則上可以對(duì)多次波成像，但需要精確的、包含正確的反射界面的速度場(chǎng)和修正的成像條件。若能提供這樣的速度場(chǎng)，RTM也就沒有必要了！二次反射波(繞射波)目前被用來進(jìn)行反射波層析反演，估計(jì)背景速度的擾動(dòng)量，提供更精確的偏移速度場(chǎng)。但是，僅僅考慮了Fresnel帶內(nèi)的二次散射波。

我們認(rèn)為在這3個(gè)層次上認(rèn)識(shí)地震數(shù)據(jù)，可以構(gòu)造出比較合適的地震波反演成像流程。事實(shí)上，構(gòu)建任何地震反演成像方法都要從疊前數(shù)據(jù)中所包含的信息能否實(shí)現(xiàn)反演目的出發(fā)。

3 參數(shù)空間及其中地震波傳播的抽象認(rèn)識(shí)

地震波在實(shí)際介質(zhì)中的激發(fā)、傳播和接收過程是地震波偏移成像和反演成像的基礎(chǔ)。用一個(gè)人為構(gòu)建的波傳播算子預(yù)測(cè)地震波傳播結(jié)果是地震波偏移成像和反演成像的核心內(nèi)容。地震波在實(shí)際介質(zhì)中的激發(fā)過程、傳播過程和波場(chǎng)的檢測(cè)過程是異常復(fù)雜的，而人為構(gòu)建的波傳播算子比較簡(jiǎn)單，僅僅描述波傳播過程，且只能描述復(fù)雜的物理波場(chǎng)的一部分，其它的波場(chǎng)成分要么處理掉，要么留下來增加偏移成像噪聲或者增加反演成像的不唯一性。這是地震勘探中地震數(shù)據(jù)處理流程復(fù)雜的基本原因。

我們把地震勘探所面對(duì)的地下介質(zhì)抽象為橫向緩慢變化的廣大沉積層中分布著由于火山活動(dòng)、后期地質(zhì)構(gòu)造運(yùn)動(dòng)和其它地球動(dòng)力運(yùn)動(dòng)所產(chǎn)生的小尺度速度異常體，譬如地震勘探中的繞射地質(zhì)體。像斷層、裂縫、地層尖滅、粗糙界面、孔洞等(大)小尺度的地質(zhì)異常體是常見的油氣運(yùn)移通道和/或儲(chǔ)集體，對(duì)它們的刻畫和描述是油氣勘探的重要目的。

針對(duì)這樣的介質(zhì)情況，我們提出如下的抽象認(rèn)識(shí)：地下介質(zhì)的速度場(chǎng)可以分為背景速度加反射界面處的速度擾動(dòng)?？捎萌缦鹿奖硎?

(3)

式中：v(x,y,z)是全波數(shù)帶的速度場(chǎng);vB(x,y,z)是其中的背景部分，δvB(x,y,z)是背景速度的擾動(dòng)量;R(x,y,z;γ,φ)是方位張角反射系數(shù)。

地震數(shù)據(jù)中不同偏移距接收到的反射同相軸的到達(dá)時(shí)間取決于背景速度vB(x,y,z)或vB(x,y,z)+δvB(x,y,z);同相軸上反射波的振幅取決于反射界面處波阻抗的變化R(x,y,z;γ,φ)。不同偏移距的數(shù)據(jù)包括透射波數(shù)據(jù)(直達(dá)波和潛波)和反射波數(shù)據(jù)(一次反射波和多次反射波)。上述觀察決定了地震速度反演的基本方法。

地震波在這樣的介質(zhì)中的傳播至少可以用兩種模式來描述：WRW模式和Born近似模式。WRW模式借用了Berkhout[19]常用的描述波傳播的方法，它認(rèn)為地下介質(zhì)是光滑的背景介質(zhì)加上具有一定反射系數(shù)的反射界面組成的。波在這樣的介質(zhì)中傳播，地面觀測(cè)到的數(shù)據(jù)可以由WRW模式對(duì)應(yīng)的正問題來描述，即下行波場(chǎng)遇到反射界面反傳回地表形成地面記錄。原則上，這樣的模式還可以描述多次波傳播。Born近似模型視地下介質(zhì)為背景速度加上小尺度的散射體，地表數(shù)據(jù)是所有地下散射體引起的散射波疊加形成的，可以僅僅包括一次散射，也可以同時(shí)包括多次散射。很顯然，WRW模式更符合目前反射地震勘探的實(shí)際，這也是很多人把勘探地震學(xué)稱為反射地震學(xué)的原因。但是，根據(jù)前面對(duì)實(shí)際地質(zhì)介質(zhì)的抽象，地下介質(zhì)中還包括很多孤立的散射點(diǎn)(體)，WRW模型需要修正以便描述反射和繞射同時(shí)存在的情況下介質(zhì)參數(shù)變化和對(duì)應(yīng)的波場(chǎng)變化之間的關(guān)系。原則上，Born近似模型包含了WRW模型，但是，Born近似條件在實(shí)際應(yīng)用中很難滿足。本質(zhì)上講，我們還需要更好的模型來描述地下介質(zhì)參數(shù)的變化和對(duì)應(yīng)的波場(chǎng)變化之間的關(guān)系。

無論如何，當(dāng)前油氣地震勘探的理論主要還是建立在反射波的基礎(chǔ)上。CMP道集處理、角度成像道集處理、AVA分析都隱含了反射波成像處理的假設(shè)，也隱含了地震波高頻近似成立的假設(shè)。

基于上述參數(shù)分解方式及對(duì)應(yīng)的波傳播模式，可以構(gòu)建出合理的地震波反演成像的實(shí)用化流程。

4 地震波反演成像的概率論觀點(diǎn)及相應(yīng)解法

綜合各種文獻(xiàn)和我們的研究經(jīng)驗(yàn)，認(rèn)為地震波反演成像技術(shù)建立在Bayes估計(jì)的理論框架上時(shí)，其理論基礎(chǔ)最為扎實(shí)，對(duì)反演問題的認(rèn)識(shí)和理解最為透徹。

首先，觀測(cè)到的地震數(shù)據(jù)是隨機(jī)信號(hào)，它是滿足一定的概率分布的。實(shí)際上，假設(shè)地震波在實(shí)際介質(zhì)中的傳播的確可以由一個(gè)正問題比較精確地預(yù)測(cè)，預(yù)測(cè)結(jié)果與實(shí)測(cè)數(shù)據(jù)的殘差也是隨機(jī)信號(hào)。地震信號(hào)的隨機(jī)性可以理解為確定性的傳播波場(chǎng)加上不可預(yù)測(cè)的隨機(jī)噪聲。一般地，我們認(rèn)為這部分不可預(yù)測(cè)的隨機(jī)噪聲滿足高斯分布。因此，隨機(jī)地震信號(hào)可以表示為：

(4)

式中：dobs為實(shí)測(cè)波場(chǎng);F(·)為地震波正演算子;m為要反演的介質(zhì)參數(shù);η為隨機(jī)噪聲。(4)式是一個(gè)非線性方程組。

事實(shí)上，可以用Newton法直接解(4)式進(jìn)行地震參數(shù)的反演，但是，這樣做不利于認(rèn)識(shí)地震參數(shù)反演問題的本質(zhì)。地震波反演成像問題不是一個(gè)簡(jiǎn)單的、可以歸結(jié)為非線性方程組求解的數(shù)學(xué)問題！

借助于Bayes估計(jì)的觀點(diǎn)，對(duì)地震波反演成像的思想的理解會(huì)更為本質(zhì)和深刻。Bayes估計(jì)的觀點(diǎn)為：參數(shù)m的反演結(jié)果或估計(jì)結(jié)果應(yīng)該由后驗(yàn)概率密度函數(shù)P(m|dobs)所定義的均值決定，即

(5)

條件概率密度P(m|dobs)可以理解為：觀測(cè)數(shù)據(jù)為實(shí)測(cè)數(shù)據(jù)dobs時(shí)對(duì)應(yīng)不同的參數(shù)m的概率?？梢?，用(5)式進(jìn)行參數(shù)m的估計(jì)是非常合乎邏輯的。但是，高維參數(shù)反演情況下的后驗(yàn)概率密度函數(shù)P(m|dobs)幾乎是無法實(shí)際計(jì)算的。

首先，利用著名的Bayes公式把后驗(yàn)概率密度函數(shù)P(m|dobs)的計(jì)算轉(zhuǎn)化為先驗(yàn)概率密度函數(shù)的計(jì)算。即:

(6)

(6)式中引入了實(shí)測(cè)數(shù)據(jù)dobs對(duì)反演參數(shù)m的先驗(yàn)概率密度函數(shù)P(dobs|m),反演參數(shù)m本身的先驗(yàn)概率密度函數(shù)P(m)和實(shí)測(cè)數(shù)據(jù)dobs的先驗(yàn)概率密度函數(shù)P(dobs)。實(shí)測(cè)數(shù)據(jù)dobs對(duì)反演參數(shù)m的先驗(yàn)概率密度函數(shù)P(dobs|m)反映了用任一個(gè)反演參數(shù)m借助正演預(yù)測(cè)算子預(yù)測(cè)實(shí)測(cè)數(shù)據(jù)時(shí)預(yù)測(cè)正確的概率或預(yù)測(cè)誤差滿足的概率。反演參數(shù)m的先驗(yàn)概率密度函數(shù)P(m)的引入相當(dāng)于引入了對(duì)反演參數(shù)的某種正則化約束。關(guān)于實(shí)測(cè)數(shù)據(jù)dobs的先驗(yàn)概率分布P(dobs)，我們一般認(rèn)為是沒有任何先驗(yàn)信息的，因此假定該先驗(yàn)概率分布是均勻分布，取常數(shù)值。

從統(tǒng)計(jì)意義上看,后驗(yàn)概率密度函數(shù)P(m|dobs)取最大值時(shí)對(duì)應(yīng)的參數(shù)m就是所要的反演結(jié)果也是很有道理的。但是，由于(6)式中概率密度函數(shù)的形式復(fù)雜，直接對(duì)(6)式進(jìn)行優(yōu)化求解是非常不方便的。

因此，一般假設(shè)P(dobs|m)和P(m)這兩個(gè)概率密度函數(shù)都是高斯函數(shù)，結(jié)合對(duì)P(dobs)的假設(shè)，有:

(7)

其中,C為常數(shù),S(m)有如下的形式:

(8)

S(m)被稱為代價(jià)函數(shù)。式中：CD是數(shù)據(jù)自相關(guān)矩陣;CM是模型自相關(guān)矩陣?？梢悦黠@看出,(8)式與Tikhonov正則化反演公式非常類似。說明Bayes估計(jì)理論具有更高的理論包容性，是反演理論更為完整的框架。由于指數(shù)函數(shù)的單調(diào)性質(zhì)，當(dāng)S(m)達(dá)到極小時(shí)，后驗(yàn)概率密度函數(shù)P(m|dobs)取最大值。因此，反演問題轉(zhuǎn)化為求如下優(yōu)化問題:

(9)

其中,H1代表參數(shù)函數(shù)取值空間,H2代表數(shù)據(jù)函數(shù)取值空間。

給出(9)式的具體求解方法之前，先分析實(shí)際地震數(shù)據(jù)反演與上述Bayes反演理論假設(shè)之間的差距。

Bayes反演理論假設(shè)實(shí)測(cè)地震波場(chǎng)是可以被預(yù)測(cè)的，預(yù)測(cè)結(jié)果與實(shí)測(cè)數(shù)據(jù)的殘差滿足高斯分布，也可以近似滿足Gauss分布，即滿足一類廣義Gauss分布。但是，實(shí)際情況下，地震波正演模擬不可能很精確地模擬實(shí)際的地震波場(chǎng)。殘差中既有觀測(cè)隨機(jī)噪聲(η)，又有不能模擬的波場(chǎng)成分，后者有可能大于前者。這會(huì)導(dǎo)致殘差滿足或近似滿足Gauss分布的假設(shè)失效，反演結(jié)果可能嚴(yán)重偏離實(shí)際情形。這是我們提出用特征波場(chǎng)或當(dāng)前的地震波正演方法能較精確模擬的波場(chǎng)進(jìn)行地震波反演的根本原因。地震子波的未知和空變也可以歸結(jié)為地震波正演方法不能很好地模擬實(shí)測(cè)波場(chǎng)，但解決的方法是盡量不用振幅逼近，而是用波場(chǎng)相關(guān)或旅行時(shí)逼近來回避子波未知導(dǎo)致的正演波場(chǎng)振幅與實(shí)測(cè)波場(chǎng)振幅的不匹配。另一方面，實(shí)測(cè)波場(chǎng)的振幅影響因素非常復(fù)雜，涉及到震源與介質(zhì)的耦合、實(shí)際震源的方向特征、檢波器與介質(zhì)的耦合、檢波器本身的響應(yīng)特征以及波在復(fù)雜介質(zhì)中的傳播過程等，并非解波動(dòng)方程可以預(yù)測(cè)的。這也是理論模型數(shù)據(jù)上經(jīng)典的全波形反演做得很好，而實(shí)際數(shù)據(jù)總做不好的原因。

地震波反演成像包含4個(gè)關(guān)鍵環(huán)節(jié)：地震波場(chǎng)數(shù)值模擬、模型表達(dá)(或模型參數(shù)化)、數(shù)據(jù)表達(dá)(或數(shù)據(jù)特征化)、反演方法。地震波數(shù)值模擬方法及線性化方法是反演過程中的關(guān)鍵問題。

(9)式描述的最優(yōu)化問題的求解，從根本上講是通過正問題的線性化把它變成一個(gè)真正的近似凸的或凸的優(yōu)化問題。這應(yīng)該是反問題最本質(zhì)的事情[20]。一旦變成了近似凸的或凸的優(yōu)化問題，其求解方法總是存在的。

將(9)式轉(zhuǎn)化為近似凸的或凸的優(yōu)化問題并沒有現(xiàn)成的方法，我們認(rèn)為從低尺度逐步過渡到高尺度的反演方法、利用特征波場(chǎng)的方法、把炮集中的波場(chǎng)分為透射波、反射波串聯(lián)反演的方法都是化非凸優(yōu)化為近似凸優(yōu)化的思路和方法。

優(yōu)化問題非凸，即存在多解性的根源，在于參數(shù)變化與波場(chǎng)變化之間的非線性性。該非線性性有兩個(gè)層次，一個(gè)是實(shí)際物理層面上的，另一個(gè)是建立的數(shù)學(xué)物理模型層次上的。一般地，我們也不去區(qū)分這兩個(gè)層次的差別。對(duì)非線性性的數(shù)學(xué)物理模型的線性化，就認(rèn)為把一個(gè)非線性的反演問題化成了線性的，非凸優(yōu)化問題就可以化為凸優(yōu)化問題。實(shí)際上，我們構(gòu)建的最小二乘反演問題、層析反演問題看起來都是凸優(yōu)化問題，但其多解性依然很強(qiáng)，主要原因是數(shù)學(xué)模型層次上的線性化并不等價(jià)于實(shí)際物理層次上的線性化。

一般地，我們只解近似凸優(yōu)化問題。若S(m)是凸的，最速下降方法、共軛梯度方法、Gauss-Newton方法和全Newton方法都可以較好地實(shí)現(xiàn)反演方法的求解。如果S(m)不是嚴(yán)格凸的函數(shù)，則S(m)的極小值點(diǎn)不再唯一，解決此類問題可以用信賴域方法、Monte Carlo方法[21]等。這些方法都可以用來求解弱非線性最優(yōu)化問題，但由于Monte Carlo方法計(jì)算量大，如果沒有一個(gè)比較好的策略很難實(shí)現(xiàn)，因此信賴域方法在這種情況下是一個(gè)比較好的選擇。

對(duì)應(yīng)非線性性很強(qiáng)的問題，難以找到一個(gè)合適的局部點(diǎn)(初始點(diǎn))把正演問題F(m)進(jìn)行線性展開得到凸性較好的S(m)。此類反演問題目前是無法求解的[22-24]。地震勘探中，小尺度異常體發(fā)育的、復(fù)雜構(gòu)造探區(qū)的反演成像連構(gòu)造都無法確定就是此類問題。針對(duì)此類問題，能做的工作依然是不斷增加數(shù)據(jù)和參數(shù)兩方面的、更豐富的信息。

對(duì)反演問題最后結(jié)果的評(píng)價(jià)原則上應(yīng)該是統(tǒng)計(jì)決策問題，主要評(píng)價(jià)依據(jù)是反演結(jié)果的均值和方差，其中方差的意義尤其重要，它反映了反演結(jié)果的可靠程度。

事實(shí)上，Bayes參數(shù)估計(jì)是一套概念，不是一種技術(shù)。在此概念下，地震波反演成像結(jié)果的好壞由很多因素決定。對(duì)Bayes框架下地震波反演成像的本質(zhì)問題的深入理解并靈活運(yùn)用才是關(guān)鍵所在。

5 FWI反演成像實(shí)用化流程——逐級(jí)線性化反演

從前面的討論可知，經(jīng)典的FWI理論是精致而完善的，具體實(shí)現(xiàn)也不復(fù)雜。但是，F(xiàn)WI反演結(jié)果總不如理論預(yù)期。主要的問題前面已經(jīng)述及，即疊前數(shù)據(jù)的不完備、地震波正演模擬方法的不完善、初始模型不夠精確、地震子波的未知和空變。因此，目前關(guān)于FWI的大量研究都是試圖讓FWI技術(shù)能在油氣勘探的實(shí)踐中發(fā)揮作用。

地震波反演成像的困難從根本上講是觀測(cè)數(shù)據(jù)與反演參數(shù)之間的線性性不強(qiáng)，以及參數(shù)的擾動(dòng)不能引起觀測(cè)數(shù)據(jù)的變化(Frechet微商等于或接近于0)。前者說明反演問題解的非唯一性，后者說明地震反演成像并非對(duì)任何參數(shù)都可以進(jìn)行。目前FWI技術(shù)發(fā)展側(cè)重的幾個(gè)方面是：①修改數(shù)據(jù)殘差和模型差的先驗(yàn)概率分布，不再要求嚴(yán)格滿足Gauss分布;②不用數(shù)據(jù)差而是用數(shù)據(jù)相關(guān)定義泛函，強(qiáng)化旅行時(shí)影響，弱化子波影響和振幅匹配;③不用原反射子波，改用反射子波的包絡(luò)來定義逼近泛函，弱化對(duì)數(shù)據(jù)中低頻的依賴性，弱化對(duì)初始模型精確性的依賴;④不用全波場(chǎng)逼近，分波型進(jìn)行FWI反演，甚至退化為分波型的基于包絡(luò)、相位、旅行時(shí)的反演。

經(jīng)典的FWI反演技術(shù)退化到分波型反演時(shí)，即認(rèn)為FWI反演可以按一個(gè)流程來實(shí)現(xiàn)。原則上，(6)式定義的經(jīng)典的FWI反演是完全自動(dòng)化的反演成像技術(shù)，不需要對(duì)疊前數(shù)據(jù)進(jìn)行任何預(yù)處理。實(shí)際上，對(duì)于復(fù)雜的疊前數(shù)據(jù)情況和復(fù)雜的地表情況，F(xiàn)WI技術(shù)要應(yīng)用于實(shí)際數(shù)據(jù)，必須構(gòu)建合適的反演成像處理流程。

根據(jù)我國地震數(shù)據(jù)的實(shí)際情況，我們提出如下實(shí)用化的流程(稱其為“FWI反演成像實(shí)用化流程”)。

1) 利用透射波(折射波和潛波)估計(jì)背景速度vB(x,y,z)。

利用中、近偏移距數(shù)據(jù)的初至波旅行時(shí)層析可以得到近地表速度，利用中、遠(yuǎn)偏移距的低頻波場(chǎng)(早至波波場(chǎng))可以估計(jì)中、深層背景速度。

2) 利用反射波和vB(x,y,z)確定反射系數(shù)R(x,y,z;γ,φ)。

利用一次(或多次)反射和散射波進(jìn)行最小二乘偏移估計(jì)反射系數(shù)，其中包括積分法的最小二乘Kirchhoff疊前深度偏移、最小二乘高斯束疊前深度偏移和波動(dòng)方程類的最小二乘波動(dòng)方程疊前深度偏移、最小二乘逆時(shí)偏移。

3) 利用反射波和vB(x,y,z)，R(x,y,z;γ,φ)，確定δvB(x,y,z)，并更新背景速度，即vB(x,y,z)?vB(x,y,z)+δvB(x,y,z)。

利用背景速度和步驟2)估計(jì)的反射系數(shù)反偏移產(chǎn)生一次反射波。通過與實(shí)測(cè)的一次反射波逼近進(jìn)行波動(dòng)理論層析，也可稱反射波FWI，估計(jì)背景速度的擾動(dòng)δvB(x,y,z)。此時(shí)同樣可以利用全波形(相位或到達(dá)時(shí))或包絡(luò)波形(包絡(luò)相位)等。

其中，步驟2)和步驟3)需要迭代執(zhí)行，直到方位張角成像道集整體拉平為止。在此基礎(chǔ)上，利用小角度反射系數(shù)可以進(jìn)入波阻抗(或密度)反演階段。

上述3步流程中的核心技術(shù)都是一種線性化的反演方法，我們把這樣的FWI實(shí)現(xiàn)方式稱為逐級(jí)線性化反演或串聯(lián)線性化反演。一般地，每項(xiàng)技術(shù)都可以進(jìn)一步利用特征波的概念來實(shí)現(xiàn)，把復(fù)雜的波進(jìn)一步分為單個(gè)震相，增加它們的實(shí)用化程度。

與上述流程相對(duì)應(yīng)，還存在一個(gè)目前石油工業(yè)界正廣泛使用的成像處理流程，產(chǎn)生保真的方位角度成像道集是其核心。該流程稱為“常規(guī)疊前深度偏移成像處理流程”，可以簡(jiǎn)述如下：

2) 水平層狀介質(zhì)假設(shè)下動(dòng)校正與速度分析技術(shù)——估計(jì)背景速度vB(x,y,z);

3) 疊前偏移剩余速度分析——估計(jì)速度擾動(dòng)量δvB(x,y,z);

4) 角度成像道集層析速度反演——估計(jì)速度擾動(dòng)量δvB(x,y,z);

5) 波動(dòng)方程疊前深度偏移——估計(jì)方位角度反射系數(shù)R(x,y,z;γ,φ)。

值得注意的是，上述兩個(gè)流程都可以包容各向異性參數(shù)、甚至Q值的估計(jì)。

對(duì)比上述兩個(gè)流程：二者的目標(biāo)是一致的，都是希望得到全(寬)波數(shù)帶的速度成分;二者的重要差異在于對(duì)背景速度的估計(jì)上。FWI反演成像實(shí)用化流程用透射波(直達(dá)波、折射波和潛波)估計(jì)背景速度。存在低頻長(zhǎng)偏移距觀測(cè)數(shù)據(jù)時(shí)，炮集中早至波成分主要是折射波和潛波，它們來自中深層。利用這些波現(xiàn)象進(jìn)行FWI可以估計(jì)背景速度vB(x,y,z)，且所建立的背景速度的精度要高于常規(guī)疊前深度偏移成像處理流程。用低頻長(zhǎng)偏移距的透射波估計(jì)背景速度是當(dāng)前FWI對(duì)實(shí)際數(shù)據(jù)處理的切實(shí)貢獻(xiàn)。在很多海上數(shù)據(jù)的FWI處理中，都是因?yàn)楸尘八俣染鹊奶岣吒纳屏顺上褓|(zhì)量(但改善并不顯著)。這就容易令人產(chǎn)生錯(cuò)覺：FWI是為了提高偏移速度精度而存在的。而實(shí)質(zhì)上，F(xiàn)WI的目標(biāo)是得到全(寬)波數(shù)帶的參數(shù)估計(jì)，它試圖提供高精度的彈性參數(shù)估計(jì)結(jié)果，直接進(jìn)入巖性儲(chǔ)層的解釋。當(dāng)然，當(dāng)前的常規(guī)疊前深度偏移成像處理流程本質(zhì)目的也是得到全(寬)波數(shù)帶的參數(shù)估計(jì)，但它往往定位在提供保真的方位角度成像道集上，基于此進(jìn)入巖性儲(chǔ)層的解釋階段。

6 陸上地震數(shù)據(jù)反演成像技術(shù)發(fā)展的關(guān)鍵問題

針對(duì)我國陸上地震數(shù)據(jù)的特點(diǎn)，推進(jìn)反演成像技術(shù)發(fā)展的關(guān)鍵在于以下幾點(diǎn)。

1) 減弱地表因素空變引起的子波特征空變。

陸上地表?xiàng)l件橫向變化劇烈，導(dǎo)致激發(fā)接收因素空間變化，反映到每一個(gè)地震道上，接收到的反射子波受炮檢點(diǎn)地表因素的影響。

減弱空變子波特征的影響首先是在估計(jì)背景速度的變化時(shí)充分利用到達(dá)時(shí)信息和相位信息，少用或不用振幅信息。其次是做好子波空變特征的一致性處理。

2) 充分利用初至波和早至波信息。

陸上數(shù)據(jù)的反射波和繞射波信噪比低，初至波和早至波一般而言具有更高的信噪比。但是，目前往往僅僅有中、小偏移距的初至波到達(dá)時(shí)被用來進(jìn)行近地表速度層析，而具有較高信噪比的初至波和早至波數(shù)據(jù)沒有得到充分的研究和利用。

3) 背景速度估計(jì)與建模。

陸上數(shù)據(jù)一般缺乏低頻長(zhǎng)偏移距數(shù)據(jù)，很難期望用透射波FWI估計(jì)背景速度vB(x,y,z)。僅僅靠動(dòng)校正和均方根速度估計(jì)在信噪比低和橫向速度變化大的情況下估計(jì)背景速度，其精度是遠(yuǎn)遠(yuǎn)不夠的，這會(huì)導(dǎo)致后續(xù)的背景速度擾動(dòng)δvB(x,y,z)估計(jì)不收斂，無法提高成像質(zhì)量。

背景速度vB(x,y,z)的估計(jì)是一個(gè)非線性反演過程。我們提出用基于Kirchhoff PSDM和成像道集拉平的Monte Carlo反演方法估計(jì)橫向變速情況下的背景速度。在2014年SEG年會(huì)中，Sajeva等用Monte Carlo方法反演背景速度為后續(xù)的FWI提供初始速度模型[25]。

4) 做好角度成像道集的層析速度反演和建模。

角度成像道集中射線層析與建模是當(dāng)前深度速度建模的主力工具。但是，它的估計(jì)精度與正則化思想與方法是密切相關(guān)的。最近幾年因?yàn)楝F(xiàn)代圖像處理技術(shù)的引入，反射界面有機(jī)地融入反演速度場(chǎng)中，極大地提高了估計(jì)精度。

對(duì)于陸上地震數(shù)據(jù)的角度成像道集的層析速度反演，在提高角度道集質(zhì)量、自動(dòng)拾取剩余時(shí)差、自動(dòng)反射界面解釋及估計(jì)傾角、Beam層析、正則化和預(yù)條件等方面的改進(jìn)會(huì)顯著提高估計(jì)精度。另外，角度成像道集的層析速度反演與各向異性參數(shù)和Q值反演可以方便地融為一體，實(shí)現(xiàn)多參數(shù)聯(lián)合建模。

5) 做好最小二乘反演成像。

在上述幾步工作的基礎(chǔ)上，疊前深度偏移成像的目標(biāo)就是至少做到生成保真的方位角度反射系數(shù)道集。我們認(rèn)為，從實(shí)用的角度看，Kirchhoff積分角度域最小二乘偏移成像和單向波最小二乘偏移成像是比較具有廣泛實(shí)用意義的。最小二乘RTM的計(jì)算量太大，難以大規(guī)模地應(yīng)用于實(shí)際數(shù)據(jù)的成像處理中。

6) 做好波阻抗反演工作。

波阻抗是與巖性變化關(guān)系更密切的參數(shù)。地震波反演成像的目標(biāo)至少應(yīng)提高到對(duì)波阻抗的反演成像上。簡(jiǎn)單地看，波阻抗與反射系數(shù)是積分和微分的關(guān)系。后者是前者的微分，反之，前者是后者的積分。從FWI反演成像角度看，小角度的反射系數(shù)與波阻抗的關(guān)系更為密切。Zhang等[26]和Duan等[27]在SEG論文中提及用小角度反射系數(shù)估計(jì)波阻抗的方法。

上述6個(gè)關(guān)鍵點(diǎn)是針對(duì)陸上地震數(shù)據(jù)進(jìn)行反演成像提出的，也可以把它們組成一個(gè)流程。與FWI反演成像實(shí)用化流程相比，它更有可能促進(jìn)反演成像方法技術(shù)的普及應(yīng)用。

另外，特別值得關(guān)注的是，陸上地震數(shù)據(jù)噪聲遠(yuǎn)大于海上地震數(shù)據(jù)，波場(chǎng)特征分解的思想更應(yīng)體現(xiàn)在陸上地震數(shù)據(jù)的FWI反演成像過程中，一是提高信噪比，二是降低正演模擬復(fù)雜波場(chǎng)的難度(特征波場(chǎng)正演模擬比全波場(chǎng)模擬簡(jiǎn)單)。另一方面，背景部分的速度估計(jì)并不需要像疊前深度偏移那樣要求數(shù)據(jù)很規(guī)則，利用波場(chǎng)的特征波分解可以把高信噪比的部分?jǐn)?shù)據(jù)分出來，僅僅用這部分?jǐn)?shù)據(jù)更新背景速度的擾動(dòng)。

信號(hào)的特征提取對(duì)整個(gè)信號(hào)處理的諸多方面(譬如數(shù)據(jù)規(guī)則化、數(shù)據(jù)壓縮、數(shù)據(jù)采樣、去噪等)都有重大影響，波場(chǎng)的特征分解在反演成像中同樣會(huì)非常重要。2014年SEG年會(huì)中已有作者關(guān)注波場(chǎng)特征分解與FWI反演數(shù)據(jù)預(yù)處理的關(guān)系[28]。

陸上地震數(shù)據(jù)的FWI反演成像需要構(gòu)建更現(xiàn)實(shí)的技術(shù)流程。若上述6個(gè)問題解決得比較好，而且每個(gè)階段的反演方法中正則化思想運(yùn)用得當(dāng)，陸上地震數(shù)據(jù)的FWI反演成像會(huì)逐步向前推進(jìn)。

7 地震波反演成像中的正則化與預(yù)條件

地震波反演成像的不唯一性和不穩(wěn)定性是其典型的特征，得到穩(wěn)定解的基本方法是正則化。正則化是地震波反演成像實(shí)用化過程中必須考慮的問題。

理論上，地震波反演成像只能進(jìn)行線性反演或弱非線性反演(可以線性化的反演)。不失一般性，可以對(duì)非線性的波場(chǎng)和參數(shù)之間的關(guān)系進(jìn)行如下線性化處理:

(10)

(10)式可以簡(jiǎn)寫為矩陣方程：

Ax=b

(11)

(11)式用來估計(jì)在m0點(diǎn)附近的δm。其中，x=δm，b=F(m)-F(m0)是數(shù)據(jù)殘差。

由于數(shù)據(jù)和稀疏矩陣中都存在誤差，(11)式本身具有欠定或超定性質(zhì)，矩陣A不滿秩,一般地，我們把對(duì)(11)式的求解化為二次型泛函的求解問題：

(12)

其中，R是正則化算子。從前面反演問題的概率論認(rèn)識(shí)中可以看出，對(duì)于(8)式右端項(xiàng)中的數(shù)據(jù)誤差項(xiàng)和模型誤差項(xiàng)，都希望是滿足高斯分布的。我們知道，當(dāng)正演算子能模擬實(shí)測(cè)信號(hào)時(shí)，數(shù)據(jù)項(xiàng)是比較符合高斯分布的，但模型項(xiàng)不容易符合高斯分布。為此，我們需要對(duì)模型進(jìn)行預(yù)處理，基本的思想是對(duì)它施加粗化算子，使它接近符合高斯分布，Laplace算子就是其中一種選擇。

從另一個(gè)角度看，反問題求解時(shí)，參數(shù)的光滑部分(低波數(shù)部分)是容易首先求解的。那么我們可以引入一個(gè)光滑算子S，作用于原參數(shù)上，即：

y=Sx

(13)

y是參數(shù)的光滑部分。把(13)式代入(11)式得到:

AS-1y≈b

(14)

可以看出估計(jì)y要容易得多，非唯一性也降低了很多。

另一方面，不同數(shù)據(jù)分量的可靠性不同，需要引入加權(quán)系數(shù)。假設(shè)加權(quán)矩陣算子為L(zhǎng)，代入(14)式得到：

LAS-1y≈Lb

(15)

(15)式實(shí)質(zhì)上是引入了對(duì)于數(shù)據(jù)的正則化處理。

對(duì)比(8)式可以看出，關(guān)于模型的正則化和數(shù)據(jù)的正則化，實(shí)質(zhì)上是試圖彌補(bǔ)(12)式中丟棄了模型協(xié)方差陣和數(shù)據(jù)協(xié)方差陣的影響。對(duì)協(xié)方差陣的地球物理意義的理解可以幫助更好地選擇(13)式和(15)式中的正則化算子。

假如對(duì)矩陣A進(jìn)行SVD分解(A=U∑VT)，并代入(15)式，可以得到：

LU∑VTS-1y=Lb

(16)

其中，U為左奇異矩陣，V為右奇異矩陣。

為了對(duì)比，(11)式引入奇異值分解結(jié)果：

U∑VTx=b

∑VTx=UTb

(17)

對(duì)比(16)式和(17)式，可知正則化算子L和S-1的選擇更基本的含義是對(duì)數(shù)據(jù)和參數(shù)引入不同的變換矩陣，使得數(shù)據(jù)中的特征成分和參數(shù)中的特征成分建立起一對(duì)一的、線性的關(guān)系。

這是最徹底的數(shù)據(jù)特征化表示和模型特征化表示。盡管上述思想不能用在大規(guī)模地震反演成像中，但它的啟發(fā)意義是深遠(yuǎn)的。對(duì)數(shù)據(jù)、對(duì)參數(shù)進(jìn)行多尺度特征分解，在此基礎(chǔ)上進(jìn)行反演是提高反演過程穩(wěn)定性的重要思想和手段。

針對(duì)(12)式，梯度導(dǎo)引類的迭代反演方法的預(yù)條件主要是對(duì)梯度場(chǎng)進(jìn)行的，其中包括施加(近似)Hessian矩陣逆的方法，也有直接對(duì)梯度進(jìn)行各種預(yù)處理的方法，這里不再贅述。TV正則化方法是十分值得關(guān)注的，很多文獻(xiàn)上有具體的討論，此處也不再贅述[23]。

8 結(jié)束語

地震波反演成像的基礎(chǔ)是多次覆蓋的數(shù)據(jù)體。目前疊前數(shù)據(jù)采集的最高技術(shù)是同時(shí)源激發(fā)進(jìn)行寬頻帶、寬方位和高密度的疊前數(shù)據(jù)采集。利用這樣的數(shù)據(jù)體對(duì)不特別復(fù)雜的地質(zhì)構(gòu)造進(jìn)行偏移成像可以給出高精度的偏移成像結(jié)果;但是若地質(zhì)情況特別復(fù)雜，即使較為可靠的偏移成像結(jié)果也難以得到。這主要不是偏移成像技術(shù)存在問題，而是估計(jì)精細(xì)偏移速度的反演成像技術(shù)存在問題。

數(shù)據(jù)體包含的基本信息是時(shí)距關(guān)系決定的波前面的到達(dá)時(shí)和波阻抗變化引起的波前面上的子波振幅變化。到達(dá)時(shí)信息是波傳播的累積效應(yīng)。波前面上的子波振幅變化是承載在到達(dá)時(shí)上的速度突變的反映。到達(dá)時(shí)對(duì)應(yīng)速度場(chǎng)的背景變化，子波振幅變化對(duì)應(yīng)速度場(chǎng)的躍變。據(jù)此，要反演的參數(shù)場(chǎng)一般可以表達(dá)為背景部分和躍變部分。對(duì)應(yīng)地，可以抽象出反射波在這樣的參數(shù)場(chǎng)中的傳播模式。

在數(shù)據(jù)表達(dá)模式、參數(shù)表達(dá)模式和波傳播模式的基礎(chǔ)上，可以很自信地構(gòu)建出適合不同數(shù)據(jù)情況的反演成像處理流程。旅行時(shí)與背景速度和背景各向異性參數(shù)之間的關(guān)系是近似線性的，可以構(gòu)建出有效的反演方法技術(shù)，一般可以滿足偏移成像方法技術(shù)的需要，給出較為精確的成像效果;地震波振幅與速度(波阻抗)擾動(dòng)之間的關(guān)系十分復(fù)雜，目前基于地震波振幅的反演結(jié)果的精度還不能滿足儲(chǔ)層描述的需要。

反演成像的目標(biāo)是得到全(寬)波數(shù)帶的參數(shù)估計(jì)結(jié)果。以逆時(shí)偏移成像(RTM)技術(shù)為標(biāo)志的、目的主要是刻畫復(fù)雜構(gòu)造界面形狀的偏移成像技術(shù)已經(jīng)發(fā)展到了接近成熟的階段。地震成像技術(shù)的發(fā)展很快會(huì)轉(zhuǎn)到以估計(jì)地下介質(zhì)彈性參數(shù)為目標(biāo)的反演成像技術(shù)上來。其中，各種層析成像技術(shù)(估計(jì)參數(shù)的背景變化)和最小二乘偏移成像技術(shù)(估計(jì)波阻抗的擾動(dòng))會(huì)成為兩項(xiàng)核心技術(shù)。

另外，把要反演的參數(shù)場(chǎng)視為一個(gè)圖像，它由紋理(高波數(shù)部分)和背景(低波數(shù)部分)組成，在反演過程中同時(shí)演化紋理的幾何和背景的大小并達(dá)到收斂，是地震速度(也包括其它參數(shù))反演的合理思想?，F(xiàn)代圖像處理中的模式提取、模式識(shí)別和模式演化會(huì)對(duì)地震波反演成像產(chǎn)生很強(qiáng)的借鑒意義和促進(jìn)作用。

[1] Tarantola A.Inversion of seismic reflection data in the acoustic approximation[J].Geophysics,1984，49(8)：1259-1266

[2] Tarantola A.Inverse problem theory and methods for model parameter estimation[M].Philadelphia,USA:SIAM,2005:1-352

[3] Kay S M.Fundamentals of statistical signal processing,volume Ⅰ:estimation theory[M].USA:Pearson Education,Inc,1993:1-625

[4] Kay S M.Fundamentals of statistical signal processing,volume Ⅱ:estimation theory[M].USA:Pearson Education,Inc,1998:1-672

[5] Beylkin G.Imaging of discontinuities in the inverse scattering problem by inversion of a causal generalized Radon transform[J].Journal of Mathematical Physics,1985,26(1)：99-108

[6] Bunks C,Salek F M,Zaleski S,et al.Multiscale seismic waveform inversion[J].Geophysics,1995,60(5)：1457-1473

[7] Lailly P.The seismic inverse problem as a sequence of before stack migrations:conference on inverse scattering[J].Expanded Abstracts of Theory and Application,Society for Industrial and Applied Mathematics,1983,206-220

[8] Mora P R.Nonlinear two-dimensional elastic inversion of multi-offset seismic data[J].Geophysics,1987,52(9)：1211-1228

[9] Mora P R.Elastic wavefield inversion of reflection and transmission data[J].Geophysics,1988，53(6)：750-759

[10] Nolet G.Seismic tomography with applications in global seismology and exploration geophysics[M].Holland:D.Reidel Publishing Company,1987:1-386

[11] Pratt R G,Worthington M H.Inverse theory applied to multisource cross-hole tomography,part I:acoustic wave-equation method[J].Geophysical Prospecting,1990,38(3):287-310

[12] Pratt R G.Inverse theory applied to multi-source cross-hole tomography,part II:elastic wave-equation method[J].Geophysical Prospecting,1990,38(3):311-330

[13] Pratt R G.Seismic waveform inversion in the frequency domain,part I:theory and verification in a physical scale model[J].Geophysics,1999,64(3):888-901

[14] Biondi B,Symes W.Angle-domain common-image gathers for migration velocity analysis by wavefield-continuation imaging[J].Geophysics,2004,69(5)：1283-1298

[15] Lambaré G,Virieux J,Madariaga R,et al.Iterative asymptotic inversion in the acoustic approximation[J].Geophysics,1992,57(9):1138-1154

[16] Shin C,Cha Y H.Waveform inversion in the Laplace-Fourier domain[J].Geophysical Journal International,2009,177(3)：1067-1079

[17] Gray S H.Seismic imaging and inversion:what are we doing,how are we doing,and where are we going?[J].Expanded Abstracts of 84thAnnual Internat SEG Mtg,2014,4416-4420

[18] Born M,Wolf E.Principles of Optics[M].7thedition.United Kingdom:Cambridge University Press,1999:1-985

[19] Berkhout A J.Combining full wavefield migration and full waveform inversion,a glance into the future of seismic imaging[J].Geophysics,2012,77(2):S43-S50

[20] Boyd S,Vandenberghe L.Convex optimization[M].United Kingdom:Cambridge University Press,2004:1-727

[21] Mosegaard K,Tarantola A.Monte Carlo sampling of solutions to inverse problems[J].Journal of Geophysical Research,1995,100(B7)：12431-12447

[22] Menke W.Geophysical data analysis:discrete inverse theory[M].USA:Academic Press,Inc,1984:1-289

[23] Vogel C R.Computational methods for inverse problems[M].Philadelphia,USA:SIAM,2002:1-183

[24] Aster R C,Borchers B,Thurber C H.Parameter estimation and inverse problems[M].USA:Elsevier Academic Press，2005:1-360

[25] Sajeva A,Bienati N,Aleardi M,et al.Estimation of velocity macro-models using stochastic full-waveform inversion[J].Expanded Abstracts of 84thAnnual Internat SEG Mtg,2014,1227-1231

[26] Zhang Y,Ratcliffe A,Duan L,et al.Velocity and impedance inversion based on true amplitude reversetime migration[J].Expanded Abstracts of 83rdAnnual Internat SEG Mtg,2013,949-953

[27] Duan L,Zhang Y,Peng L Z.Band-limited impedance perturbation inversion using cross-correlative least-squares RTM[J].Expanded Abstracts of 84thAnnu-al Internat SEG Mtg,2014,3720-3725

[28] Bi H B,Lin T.Impact of adaptive data selection on full waveform inversion[J].Expanded Abstracts of 84thAnnual Internat SEG Mtg,2014,1094-1098

(編輯：陳 杰)

Analysis of seismic inversion imaging and its technical core issues

Wang Huazhong,Feng Bo,Wang Xiongwen,Hu Jiangtao,Li Hui

(WavePhenomenaandInversionImagingResearchGroup(WPI),SchoolofOceanandEarthScience,TongjiUniversity,Shanghai200092,China)

The conventional seismic inversion imaging is divided into migration velocity analysis (MVA) and tomography,prestack migration imaging and AVA analysis/inversion.The key techniques in the three inversion imaging methods are the core technologies of modern exploration seismology.Full waveform inversion (FWI) is a complete theoretical framework of seismic inversion imaging,and which can roll up the three conventional inversion imaging methods into one and give a desired imaging result in principle.However,due to the incompleteness of pre-stack seismic data,imperfectness of forward model,inaccuracy of initial model,unknown and spatial variation of seismic wavelet,the practical problems have not yet been solved with FWI method in the strict sense.In this paper we analyze the existing problems of classical FWI method and propose a practical seismic inversion imaging workflow,that is to divide classic FWI into transmissive wave tomography,LS-PSDM imaging and reflections tomography.Aimed at the features of land seismic data,we point out that the establishment of shallow velocity model,background velocity model,velocity model for CIG tomography,LS-PSDM imaging,generation of the incident angle gather and dip angel gather,impedance inversion by small angel CIG gathers are crucial for the promotion of the application and self-development of FWI technique.

the inversion imaging theory,multistage Born approximation,full waveform inversion (FWI),stepwise linearization inversion,practical seismic inversion imaging route

2015-01-03;改回日期：2015-01-27。

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

國家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃(973計(jì)劃)(2011CB201002)、國家自然科學(xué)基金(41374117)和國家科技重大專項(xiàng)(2011ZX05005-005-008HZ，2011ZX05006-002，2011ZX05023)共同資助。

P631

A

1000-1441(2015)02-0115-11

10.3969/j.issn.1000-1441.2015.02.001