国产日韩欧美一区二区三区三州_亚洲少妇熟女av_久久久久亚洲av国产精品_波多野结衣网站一区二区_亚洲欧美色片在线91_国产亚洲精品精品国产优播av_日本一区二区三区波多野结衣 _久久国产av不卡

?

針對巖性儲(chǔ)層的定量地震波成像

2023-02-09 12:42:26王華忠吳成梁許榮偉張如一
石油物探 2023年1期
關(guān)鍵詞:子波波阻抗波數(shù)

王華忠,吳成梁,盛 燊,許榮偉,雷 霆,張如一,2

(1.波現(xiàn)象與智能反演成像研究組(WPI),同濟(jì)大學(xué)海洋與地球科學(xué)學(xué)院,上海200092;2.中石化石油物探技術(shù)研究院有限公司,江蘇南京211103)

油氣地震勘探最終目的是精確地描述油氣藏并對其進(jìn)行準(zhǔn)確的含油氣性分析,做出最佳的鉆井決策,得到最高的油氣勘探效益。油氣地震勘探核心問題是由疊前地震數(shù)據(jù)及與待估計(jì)參數(shù)相關(guān)的先驗(yàn)信息,進(jìn)行寬波數(shù)帶的彈性參數(shù)估計(jì)(或稱廣義的高精度地震波成像),與巖石物理知識(shí)結(jié)合,進(jìn)行精確的油氣藏描述和準(zhǔn)確的含油氣性評價(jià)。寬波數(shù)帶的彈性參數(shù)估計(jì)問題是基于疊前地震數(shù)據(jù)和待估計(jì)參數(shù)相關(guān)的先驗(yàn)信息,基于地震波理論和波動(dòng)方程(及其各種簡化形式),基于Bayes參數(shù)估計(jì)理論,提出并求解一個(gè)信息不足情形下的(強(qiáng))非線性反問題(系統(tǒng)參數(shù)反演問題)。高精度地震波成像的信息源是“兩寬兩高(寬方位、寬頻帶、高密度、高信噪比)”的地震數(shù)據(jù)。而僅靠大炮數(shù)據(jù)采集的“兩寬兩高”疊前地震數(shù)據(jù)還遠(yuǎn)不能滿足高精度、定量、針對巖性儲(chǔ)層描述的高精度地震波成像的需求,還需要補(bǔ)充待估計(jì)的與地下介質(zhì)巖性變化相關(guān)的先驗(yàn)信息。可以說,油氣地震勘探的整個(gè)發(fā)展歷史,就是在地震波成像理論的指導(dǎo)下,由地震數(shù)據(jù)采集技術(shù)(也包括其它先驗(yàn)信息收集技術(shù))不斷向前發(fā)展的牽引下而向前發(fā)展的。而地震數(shù)據(jù)采集技術(shù)的進(jìn)步奠基于裝備技術(shù)的發(fā)展。這些構(gòu)成了當(dāng)前油氣地震勘探領(lǐng)域的核心技術(shù)概念。

地震波成像理論和技術(shù)將疊前地震數(shù)據(jù)及與待估計(jì)參數(shù)相關(guān)的先驗(yàn)信息映射成地下介質(zhì)彈性參數(shù)、物性參數(shù),進(jìn)而實(shí)現(xiàn)高精度的油氣藏描述及高效油氣開發(fā)。

常規(guī)的地震波成像以地下界面的帶限反射系數(shù)為成像目標(biāo),高精度的速度建模為保真和高分辨的帶限反射系數(shù)成像提供更精確的背景速度場。由于震源激發(fā)端和檢波端各種因素的不可控,地震波成像方法技術(shù)并不能完全消除地震波傳播過程中的各種影響因素,即便做到保真和高分辨的帶限反射系數(shù)的成像也非常不容易。但是,深水、深層情況下碎屑巖儲(chǔ)層、火山巖儲(chǔ)層、碳酸鹽巖儲(chǔ)層和非常規(guī)儲(chǔ)層都已成為油氣勘探的常規(guī)目標(biāo)。針對這些復(fù)雜巖性油藏地震勘探的需求,地震波成像技術(shù)由定性的、帶限反射系數(shù)為成像目標(biāo)轉(zhuǎn)為定量的、寬帶波阻抗為成像目標(biāo)成為必然選擇。

當(dāng)前,石油工業(yè)界開展巖性油藏勘探的地震波成像主要依賴3項(xiàng)關(guān)鍵技術(shù),即層析速度反演技術(shù)、保真和高分辨方位角度帶限反射系數(shù)成像方法技術(shù)、AVA分析/反演方法技術(shù)。3項(xiàng)技術(shù)是串聯(lián)的,最后一項(xiàng)技術(shù)反映了角度反射系數(shù)反演界面兩邊彈性參數(shù)的變化(ΔvP為縱波速度擾動(dòng);ΔvS為橫波速度擾動(dòng);Δρ為密度擾動(dòng))。存在的主要問題是AVA分析/反演的基本理論基礎(chǔ)是Zoeppritz方程,而Zoeppritz方程是在平面波入射到廣大的平反射界面上的假設(shè)情況下推導(dǎo)出的。在碳酸鹽巖油藏、火山巖油藏等情形下是不合適的假設(shè)。另一方面,保真帶限反射系數(shù)估計(jì)本身就很困難,反演結(jié)果存在很大的不確定性。常規(guī)1D波阻抗反演是當(dāng)前波阻抗建模的主要方法技術(shù),但是1D波阻抗反演僅適用于水平層狀介質(zhì)情形[1]。歷史上,該項(xiàng)技術(shù)主要用于緩橫向變化的碎屑巖沉積地層下的油氣勘探中,顯然已經(jīng)完全不適用于目前橫向變化劇烈的復(fù)雜油氣藏的勘探。

Bayes參數(shù)估計(jì)理論下的全波形反演(FWI)技術(shù)已經(jīng)成為業(yè)界公認(rèn)的高精度地震波成像的標(biāo)志性技術(shù)[2]。從理論本質(zhì)上看,FWI技術(shù)的目標(biāo)是估計(jì)寬波數(shù)帶的速度變化(或者波阻抗變化),它是定量的地震波成像方法。眾所周知,當(dāng)前FWI被定義為主要利用透射波估計(jì)高精度的背景速度模型,為逆時(shí)偏移(RTM)或最小二乘逆時(shí)偏移(LS_RTM)提供更好的初始速度模型。即便如此,也需要退化成利用透射波走時(shí)信息,并且在海上低頻長偏移距數(shù)據(jù)情形下,才有可能產(chǎn)生實(shí)用的效果。根本的困難在于實(shí)際觀測數(shù)據(jù)不能滿足FWI反演的理論假設(shè)以及復(fù)雜介質(zhì)變化(小尺度、強(qiáng)變化)情形下FWI反演的高度非線性性。因此,實(shí)際數(shù)據(jù)情形下利用FWI反演技術(shù)進(jìn)行寬帶波阻抗為目標(biāo)的成像顯然是不可能的。

實(shí)際數(shù)據(jù)的地震波反演成像從來不是一個(gè)純數(shù)學(xué)問題,即便“兩寬兩高”的疊前地震數(shù)據(jù)依然不能滿足巖性儲(chǔ)層定量反演成像的需求,必須發(fā)展出地質(zhì)信息約束下的地震波成像新理論與新方法。定性的地震波偏移成像技術(shù)的研究已經(jīng)趨于完善。但是,任何更為高端的技術(shù)都是在現(xiàn)有技術(shù)基礎(chǔ)上發(fā)展出來的。充分保留現(xiàn)有方法技術(shù)的優(yōu)點(diǎn),充分借鑒Bayes參數(shù)估計(jì)理論下的FWI的理論優(yōu)勢,充分利用大數(shù)據(jù)和機(jī)器學(xué)習(xí)(ML)中的理念和算法[3],發(fā)展信息融合的深度域真三維寬帶波阻抗成像,認(rèn)為是合理的針對巖性儲(chǔ)層的定量地震波成像技術(shù)路線。

寬波數(shù)帶波阻抗可以被分為低波數(shù)帶(0~3Hz)、中波數(shù)帶(4~8Hz)、高波數(shù)帶(9~50Hz)和超高波數(shù)帶(51~80Hz)4個(gè)部分[4](注:其中的數(shù)字劃分是相對的)。背景阻抗的信息來源一定是背景速度場和背景密度場;高波數(shù)帶阻抗的信息源一定是寬帶反射系數(shù)(或?qū)拵У乃俣葦_動(dòng))。中波數(shù)波阻抗的信息來源更多地依賴于與地下介質(zhì)波阻抗變化有關(guān)的(地質(zhì))先驗(yàn)信息;超高波數(shù)波阻抗的信息源應(yīng)該是測井?dāng)?shù)據(jù)和巖石物理知識(shí)。

另外,分析可知,實(shí)際數(shù)據(jù)情形下的FWI也不能得到由低波數(shù),到中波數(shù),再到高波數(shù)的波阻抗(或速度)反演結(jié)果。眾所周知的關(guān)系式FWI=Tomography+LS_RTM僅僅在反演騙局(Inverse Crime)意義下是成立的[5]。對于實(shí)際數(shù)據(jù),這個(gè)公式會(huì)引起誤導(dǎo)。在實(shí)際數(shù)據(jù)情形下,基于反射地震子波振幅的LS_RTM方法并不能得到與地下反射系數(shù)同等量級(jí)的反演結(jié)果,因此高波數(shù)速度擾動(dòng)的估計(jì)過程是不可能收斂的(連反演結(jié)果量級(jí)都不能與低波數(shù)反演結(jié)果保持一致,如何能產(chǎn)生收斂的反演結(jié)果?)。

依據(jù)上述分析,我們(WPI研究組)提出了特征波反演成像(characteristic wave inversion imaging,CWI)與基于信息融合的寬帶波阻抗建模(wide band impedance modeling,WBIM)的技術(shù)路線[6],是一種針對巖性儲(chǔ)層的定量地震波成像方法技術(shù)系列。關(guān)鍵技術(shù)包括:背景速度模型建立方法、背景密度模型建立方法、保真的寬帶反射系數(shù)估計(jì)方法和基于信息融合的寬帶波阻抗建模方法。必須說明的是,提高各環(huán)節(jié)建模精度的先驗(yàn)信息獲取方法(包括先驗(yàn)信息的正則化約束方法)也是十分重要的關(guān)鍵方法技術(shù)。原因是即便有“兩寬兩高”的疊前地震數(shù)據(jù)還是不能很好地滿足定量地震波成像的需求,必須引入盡可能多的先驗(yàn)信息,才有可能進(jìn)一步提升各波數(shù)帶的波阻抗建模精度。這與當(dāng)前大數(shù)據(jù)、人工智能和機(jī)器學(xué)習(xí)等技術(shù)越來越多地滲透到勘探地震全流程中的趨勢是完全一致的。因此,基于信息的、定量的地震波成像代表了地震波成像技術(shù)的發(fā)展趨勢。

理論數(shù)據(jù)和實(shí)際數(shù)據(jù)的應(yīng)用實(shí)踐證明了我們提出的特征波反演成像(CWI)與基于信息融合的寬帶波阻抗建模(WBIM)流程及方法技術(shù)的有效性。

1 定量地震波成像理論基礎(chǔ)

前已述及,油氣地震勘探的核心是寬波數(shù)帶的彈性參數(shù)估計(jì),它提供了高精度油藏描述的基礎(chǔ)信息,結(jié)合井?dāng)?shù)據(jù)、巖石物理和油氣地質(zhì)知識(shí),由高精度油藏描述和含流體性評價(jià)提供最佳的井位決策,獲得最佳的勘探效益。很顯然,油氣勘探中寬波數(shù)帶的彈性參數(shù)估計(jì)(或稱高精度地震波成像)占據(jù)核心地位。我們將其定位成提出并求解一個(gè)信息不足情形下的(強(qiáng))非線性反問題(系統(tǒng)參數(shù)反演問題)。

寬波數(shù)帶的彈性參數(shù)估計(jì)(或稱高精度地震波成像)的信息源包括疊前地震數(shù)據(jù)(盡可能是“兩寬兩高”的)及與待估計(jì)參數(shù)有關(guān)的先驗(yàn)信息。

Bayes估計(jì)理論奠定了寬波數(shù)帶的彈性參數(shù)估計(jì)(或稱高精度地震波成像)的理論基礎(chǔ)。疊前地震數(shù)據(jù)被記為dobs(xS,xR,t);待估計(jì)的寬波數(shù)帶的彈性參數(shù)被記為m(x)。疊前地震數(shù)據(jù)一般被認(rèn)為是五維數(shù)據(jù)體,反演成像結(jié)果是深度域真三維的。疊前地震數(shù)據(jù)被認(rèn)為是隨機(jī)過程,實(shí)際野外數(shù)據(jù)采集被認(rèn)為是該隨機(jī)過程的一次實(shí)現(xiàn)。認(rèn)為波動(dòng)方程(或某種變種)建立起了m(x)與dobs(xS,xR,t)之間一定的預(yù)測關(guān)系。很顯然,這樣的預(yù)測關(guān)系是不確定的,波動(dòng)方程(或某種變種)并不能完全地預(yù)測dobs(xS,xR,t)中包含的各種復(fù)雜的波現(xiàn)象。實(shí)測數(shù)據(jù)的不確定性,波動(dòng)方程預(yù)測結(jié)果的不確定性,使得我們有理由認(rèn)為估計(jì)的結(jié)果也是不確定的。這就構(gòu)成了Bayes估計(jì)理論的思想基礎(chǔ)。

經(jīng)典的Bayes公式為:

(1)

式中:P(m)代表先驗(yàn)的表達(dá)待估計(jì)參數(shù)的不確定性的概率密度函數(shù);P(dobs|m)代表先驗(yàn)的波動(dòng)方程(或某種變種)預(yù)測實(shí)測數(shù)據(jù)時(shí)不確定性的概率密度函數(shù)(也稱似然函數(shù)),這里面既有波動(dòng)方程的不合適引起的預(yù)測結(jié)果的不確定性,也有實(shí)際數(shù)據(jù)觀測中各種復(fù)雜因素導(dǎo)致的觀測結(jié)果的不確定性,震源激發(fā)和檢波器性能引起的不確定性也都包容在這個(gè)概率密度函數(shù)中??梢钥闯?(1)式中包含了實(shí)際地震數(shù)據(jù)觀測、波動(dòng)方程對高精度地震波成像的不滿足;也包含了引入與待估計(jì)參數(shù)有關(guān)的先驗(yàn)信息的必要性。僅僅靠波動(dòng)方程對實(shí)際數(shù)據(jù)的預(yù)測誤差很小(或似然函數(shù)P(dobs|m)取最大),并不能(遠(yuǎn)遠(yuǎn)不能)保證反演結(jié)果的可靠性。P(dobs)是描述采集數(shù)據(jù)過程表現(xiàn)出的不確定性的概率密度函數(shù)。一般地,我們很難獲取這方面的先驗(yàn)信息,索性假設(shè)P(dobs)滿足均勻分布,不再引入先驗(yàn)信息。因此,(1)被改寫為:

P(m|dobs)∝P(dobs|m)P(m)

(2)

實(shí)際上,P(m)和P(dobs|m)都是很難實(shí)際計(jì)算得到的。因此,先驗(yàn)地假定它們均滿足Gauss概率分布。很多實(shí)際隨機(jī)過程被假定滿足Gauss概率分布,這是有大數(shù)定理和中心極限定理支持的。因此,有:

(3a)

(3b)

式中:Cm代表待估計(jì)參數(shù)的協(xié)方差矩陣,M是它的維度;Cd代表觀測數(shù)據(jù)的協(xié)方差矩陣,N是它的維度。K(m)=d抽象地表示波動(dòng)方程(或其變種)對實(shí)測數(shù)據(jù)的預(yù)測關(guān)系,數(shù)學(xué)上更多用第一類Fredholm積分方程來表示,而非用微分形式的波動(dòng)方程。地震勘探中Born近似后的散射場計(jì)算式就是典型的第一類Fredholm積分方程,它建立起了地下介質(zhì)參數(shù)的擾動(dòng)量與波場擾動(dòng)量之間的數(shù)學(xué)關(guān)系式。

將(3)代入(2)式得到:

(4)

很顯然,根據(jù)指數(shù)函數(shù)的單調(diào)性,后驗(yàn)概率密度最大化等價(jià)于S(m)取極小。S(m)被稱為代價(jià)函數(shù)或誤差泛函。Bayes參數(shù)估計(jì)問題轉(zhuǎn)化為了誤差泛函求極小問題。

(5)

這就是最經(jīng)典的FWI反演理論。其目標(biāo)顯然是全頻帶的、定量的參數(shù)估計(jì)。理論上看,Bayes估計(jì)理論下的FWI似乎是最完美的地震波成像方法。

(5)式求解根本的困難在于實(shí)際觀測數(shù)據(jù)不能滿足FWI反演的理論假設(shè);在于復(fù)雜介質(zhì)變化(小尺度、強(qiáng)變化)情形下它表現(xiàn)出高度非線性性。

經(jīng)典的勘探地震中,我們并沒有硬解(5)式定義的反問題,而是把正問題(波動(dòng)方程對數(shù)據(jù)的預(yù)測關(guān)系)退化成線性關(guān)系,即

(6a)

(6b)

式中:d0=K(m0),F=?K/?m。將(6)式代入(5)式,得到:

(7)

(7)式是地震波層析成像和LS_PSDM共同遵循的參數(shù)估計(jì)反問題提法。因此有FWI=Tomography+LS_PSDM的說法。層析成像解決背景速度估計(jì)問題;LS_PSDM解決帶限反射系數(shù)(或速度擾動(dòng)量)估計(jì)問題。通過(6)式表示的線性化過程,(5)式定義的非線性反問題被轉(zhuǎn)化為線性的。線性化的反問題理論上具備凸性,存在唯一解。因此(7)式的求解理論上應(yīng)該更穩(wěn)定、迭代解的收斂性更好、反演結(jié)果更精確[7-8]。實(shí)際上,常規(guī)的地震波成像過程中的三項(xiàng)關(guān)鍵技術(shù):層析成像、偏移成像和AVA反演成像都是線性反演理論在勘探地震中的應(yīng)用實(shí)踐??梢钥闯?常規(guī)的地震波反演成像目的也是期望得到寬波數(shù)帶的彈性參數(shù)估計(jì)。但核心問題出在Zoeppritz方程的理論假設(shè)(平面波入射到廣大的地下介質(zhì)平反射界面上)不適應(yīng)當(dāng)前常見的橫向變化劇烈的介質(zhì)情形;偏移成像給出的反射系數(shù)不保真,與地下介質(zhì)反射界面的“真”反射系數(shù)差異可能很大,AVA反演缺乏可靠的數(shù)據(jù)基礎(chǔ)。常規(guī)的1D波阻抗反演本質(zhì)上是地下介質(zhì)水平層狀假設(shè)下的FWI反演,(5)式定義的FWI遇到的問題,1D FWI同樣不能回避。更重要的是,目前油氣勘探遇到的復(fù)雜介質(zhì)和復(fù)雜油藏情形已經(jīng)不符合其水平層狀介質(zhì)假設(shè)。

因此,針對巖性油藏的定量地震波成像需要新的更合理的技術(shù)路線。我們是通過修改(5)式定義的FWI形成由特征波反演成像+基于信息融合的寬帶波阻抗建模構(gòu)成的定量地震波成像路線。

由前述Bayes估計(jì)理論可知,波動(dòng)方程(或各種變種)對實(shí)測數(shù)據(jù)(中的波現(xiàn)象)的預(yù)測關(guān)系,本質(zhì)是待估計(jì)參數(shù)與波現(xiàn)象之間的關(guān)系是否接近線性決定了形成的反問題是否凸的,反問題是否有很好的凸性決定了反問題的求解是否能得到穩(wěn)定可靠的解。先驗(yàn)信息的引入可以促進(jìn)解的穩(wěn)定性和精度進(jìn)一步提升,但過分依賴初始模型或先驗(yàn)信息提高解的精度邏輯上是不合理的。因?yàn)橄闰?yàn)信息的獲取往往比“兩寬兩高”的地震數(shù)據(jù)采集更困難,代價(jià)更大。比如我們不可能鉆很多井來獲取關(guān)于地下介質(zhì)彈性參數(shù)、物性參數(shù)的先驗(yàn)信息。

為此,我們提出了特征波反演成像的概念,(5)式被改寫成:

(8)

我們認(rèn)為,特征波成像主要包括初至波層析成像和特征反射層導(dǎo)引下的特征反射波層析成像。常規(guī)的射線理論折射波層析、基于各種反射波時(shí)距關(guān)系的NMO速度分析、共成像點(diǎn)道集射線層析成像也可歸于特征波成像的范圍。

很顯然,特征波成像是基于不同波現(xiàn)象同相軸上地震子波的走時(shí)估計(jì)背景速度,這是一個(gè)定量反演過程。所得的背景速度與地下介質(zhì)速度量級(jí)是一致的。它們可被認(rèn)為是定量的成像方法。

(7)式規(guī)定的另一個(gè)線性化反演問題是最小二乘疊前深度偏移(LS_PSDM)或LS_RTM。常規(guī)的地震波偏移成像是LS_PSDM或LS_RTM的簡化實(shí)現(xiàn)。(7)式的數(shù)學(xué)求解可以表示為:

(9)

常規(guī)的地震波偏移成像的數(shù)學(xué)表述式為:

(10)

實(shí)際上,石油工業(yè)界常用的偏移成像思想和方法技術(shù)并非由上述Bayes估計(jì)理論退化產(chǎn)生的,而是由物理直覺與波傳播算子相結(jié)合導(dǎo)出。認(rèn)為地表觀測的波場是由地表主動(dòng)源產(chǎn)生的下行波激活地下介質(zhì)波阻抗異常體,產(chǎn)生上行繞射波場,所有地下介質(zhì)中的波阻抗異常體產(chǎn)生的上行波場的疊加形成地表觀測波場。尤其考慮到多次覆蓋的地表數(shù)據(jù)激發(fā)接收方式,地震波成像就是把所有地表接收的、由所有炮激發(fā)的上行繞射波場推到或收斂到產(chǎn)生它們的波阻抗異常體上。Kirchhoff積分偏移成像就是直接地按上述物理邏輯實(shí)現(xiàn)的;單平方根波動(dòng)方程偏移是按正向外推炮端激發(fā)波場與反向外推接收端波場并在波阻抗異常體處的相關(guān)(相關(guān)最大對應(yīng)偏移速度模型是正確的)進(jìn)行成像的;雙平方根波動(dòng)方程偏移是按炮點(diǎn)和檢波點(diǎn)同時(shí)沉降的方式,把炮點(diǎn)和檢波點(diǎn)外推到波阻抗異常體處按成像條件(上行波出發(fā)時(shí)刻為零)提取成像值進(jìn)行成像的。但是,這些具體的實(shí)現(xiàn)方法本質(zhì)上與Bayes估計(jì)理論退化出的(10)式是統(tǒng)一的。只不過前者物理意義和幾何意義更清晰;后者的數(shù)學(xué)算法更嚴(yán)謹(jǐn)。

很顯然,石油工業(yè)界常用的偏移成像只不過是對地下波阻抗異常體的位置進(jìn)行定位,定位精度依賴于背景速度。定位結(jié)果的表現(xiàn)是由地震子波的振幅體現(xiàn)出來的,正是“放置或搬家”到成像點(diǎn)處反射子波的振幅體現(xiàn)出地下介質(zhì)反射系數(shù)的相對變化(并期望這種相對變化能比較真實(shí)地反映地下介質(zhì)反射系數(shù)真實(shí)情況)。直到目前的油氣地震地質(zhì)解釋主要還是利用了地震波成像提供的波阻抗異常體的空間定位信息和保真的反射強(qiáng)度信息。事實(shí)上,相關(guān)成像條件得到定位結(jié)果是由下行波場與外推到成像點(diǎn)處的下行波場中子波的相關(guān)系數(shù)體現(xiàn)出來的,與地震子波的振幅關(guān)系更遠(yuǎn)。如果不強(qiáng)調(diào)偏移成像的保真處理,成像結(jié)果得到的帶限反射系數(shù)與地下介質(zhì)反射系數(shù)之間可能相差極大。

實(shí)際觀測波場中,各震源激發(fā)產(chǎn)生的地震波能量(振幅)有復(fù)雜的影響因素,炮和炮之間差異很大;檢波器測震結(jié)果的物理意義不同,有位移檢波器、速度檢波器、加速度檢波器和壓力檢波器,各檢波器與介質(zhì)的耦合相差很大,檢波器本身對不同頻率響應(yīng)的一致性很難保證,導(dǎo)致檢波器檢測的地震波振幅是相對的。而且各檢波器對同一炮激發(fā)的震動(dòng)的響應(yīng)差異很大,即便相鄰檢波器測量的同一反射層的反射子波振幅也很不相同。這些復(fù)雜的因素并不能被正問題K(m)=d所描述。地震波反演成像中的正算子很簡單,常密度聲波方程最常用,而且又被引入了Born近似,復(fù)雜介質(zhì)對波場的改造作用不能在成像過程中得到消除。因此,實(shí)際數(shù)據(jù)情況下,即便使用LS_RTM方法也得不到收斂的定量成像結(jié)果。理論框架下的FWI和LS_RTM在自身的理論假設(shè)下可以得到可靠的定量成像結(jié)果,遺憾的是,這只是個(gè)反演騙局(Inverse Crime)。

當(dāng)前,FWI_Imaging[9-11]或全波數(shù)帶FWI[12]試圖由透射波FWI估計(jì)低波數(shù)帶背景速度,然后進(jìn)到反射波FWI估計(jì)中波數(shù)帶的背景速度,再進(jìn)展到估計(jì)高波數(shù)帶的速度擾動(dòng)量(這已經(jīng)是LS_RTM解決的問題,這是LS_RTM被稱為高波數(shù)參數(shù)反射波FWI的原因)。事實(shí)上,由于FWI反演迭代過程中,由利用波現(xiàn)象同相軸上地震子波的走時(shí)信息到利用地震子波的振幅信息,迭代反演中的梯度項(xiàng)由定量變?yōu)椤岸ㄐ浴?不受控制的FWI根本不可能由定量的背景速度的更新自動(dòng)地轉(zhuǎn)到定量的高波數(shù)擾動(dòng)量的更新。因此,盡管在Inverse Crime意義下,FWI的梯度項(xiàng)由背景速度的更新自動(dòng)地轉(zhuǎn)到高波數(shù)擾動(dòng)量的更新似乎是理所當(dāng)然的(事實(shí)上,也不是那么的理所當(dāng)然);但對于實(shí)際數(shù)據(jù),這是不可能自動(dòng)發(fā)生的。

影響地震波振幅的因素非常復(fù)雜,基于地震子波振幅的參數(shù)估計(jì)或定量成像到目前為止并沒有完善的解決方案。上述分析清楚地說明,針對巖性儲(chǔ)層的定量地震波成像的真正困難在于高波數(shù)參數(shù)的定量化估計(jì)。更明確地,保真與寬帶反射系數(shù)成像是定量地震波成像的難點(diǎn)。

假定已經(jīng)獲得了定量的背景速度、背景密度、保真與寬帶反射系數(shù),基于上述3種成像結(jié)果的地震地質(zhì)解釋依然不是深度域真三維空間的油氣地質(zhì)信息的真實(shí)反映。深度域真三維空間的寬帶波阻抗才是更能反映地下介質(zhì)巖性變化情況的成像結(jié)果。深度域真三維空間的寬帶波阻抗自然不能作為FWI的成像目標(biāo),這只會(huì)徒然增加反問題的非線性性,求解過程更不收斂,多解性更嚴(yán)重。這正是我們提出信息融合寬帶波阻抗建模的本質(zhì)原因,信息融合寬帶波阻抗建模避開了FWI反演寬帶波阻抗的強(qiáng)非線性性,可以將更多的與阻抗有關(guān)的先驗(yàn)信息融入到寬帶波阻抗建模過程中,從而充分借用當(dāng)前機(jī)器學(xué)習(xí)領(lǐng)域發(fā)展的新思想和新算法。

定量地震波反演成像是工業(yè)界有急迫需求且理論上也存在待解決問題的研究方向。我們(WPI研究組)提出的CWI(Characteristic Wave Inversion)+WBIM(Wide Band Impedance Modeling)策略就是一種定量的地震成像技術(shù)路線。

2 定量地震波成像技術(shù)路線

前已述及,當(dāng)前油氣地震勘探主要還是利用了地震波成像提供的地下介質(zhì)波阻抗異常體的定位信息,以及偏移成像算法附帶產(chǎn)生的相對保真的反射強(qiáng)度信息,進(jìn)行地震地質(zhì)解釋,在井?dāng)?shù)據(jù)和油氣地質(zhì)知識(shí)的輔助下進(jìn)行井位決策和鉆井方案決策。

很顯然,地震波成像由定性的、主要反映地下介質(zhì)幾何形態(tài)附帶表現(xiàn)地下介質(zhì)反射強(qiáng)度變化的帶限反射系數(shù)成像發(fā)展到定量的、反映地下介質(zhì)巖性變化的寬帶波阻抗成像,顯然滿足了油氣工業(yè)的重大需求。圖1展示了由疊前地震數(shù)據(jù)與待估計(jì)參數(shù)相關(guān)的先驗(yàn)信息映射成油藏描述需要的定量地震成像結(jié)果的3條技術(shù)路線。

由前述理論分析知道,Bayes估計(jì)理論下的FWI是定量地震波成像最為完善的方法技術(shù)。但是,在實(shí)際數(shù)據(jù)情況下,FWI直接反演寬帶波阻抗是不可能的。這是一個(gè)高度非線性的反演問題。

常規(guī)的做法包括3項(xiàng)核心技術(shù),即背景速度層析成像;帶限的角度反射系數(shù)成像和AVA反演(也包括1D波動(dòng)方程波阻抗反演),存在的問題(分析后認(rèn)為)主要是AVA反演的Zoeppritz方程理論假設(shè)不適應(yīng)當(dāng)前復(fù)雜的勘探目標(biāo)和角度反射系數(shù)成像結(jié)果存在嚴(yán)重的不確定性。另外,1D波動(dòng)方程波阻抗反演的理論假設(shè)與當(dāng)前勘探實(shí)際差異過大,橫向變化劇烈的介質(zhì)情況下,根本不存在滿足1D波動(dòng)方程的1D數(shù)據(jù)。

可以看出,圖1中所示的走向定量油藏描述的常規(guī)技術(shù)路線越來越不適應(yīng)當(dāng)前復(fù)雜介質(zhì)與復(fù)雜油藏的勘探實(shí)際。

HAFFINGER等[13]提出基于FWI反演的AVA彈性參數(shù)反演流程(wave-equation based AVO inversion,WEB-AVO)試圖解決定量儲(chǔ)層描述問題。但是,基于FWI背景速度及疊前深度偏移成像道集和1.5D彈性參數(shù)反演這三項(xiàng)關(guān)鍵技術(shù)的WEB-AVO存在很多邏輯上令人質(zhì)疑的地方,譬如疊前深度偏移給出的零偏移距地震剖面是否滿足1.5D彈性波方程的要求?1.5D彈性波反演依然是個(gè)強(qiáng)非線性問題,如何能得到穩(wěn)定的反演解?另外,1.5D彈性參數(shù)反演假設(shè)地下介質(zhì)水平層狀也不適應(yīng)當(dāng)前橫向變化復(fù)雜的勘探目標(biāo)。

圖1 由疊前地震數(shù)據(jù)與待估計(jì)參數(shù)相關(guān)的先驗(yàn)信息映射成油藏描述需要的定量地震成像結(jié)果的3條技術(shù)路線

根據(jù)上一節(jié)的理論分析,我們(WPI研究組)提出了另外一條定量地震波成像技術(shù)路線(圖2)。

該技術(shù)路線中包括背景速度建模、背景密度建模、寬帶反射系數(shù)(或高波數(shù)參數(shù)擾動(dòng)量)估計(jì)及信息融合寬帶波阻抗建模四大部分。

背景速度和背景密度構(gòu)成的背景阻抗對于巖性油藏的描述與評價(jià)具有重要意義。背景速度模型的構(gòu)建包含如下5項(xiàng)關(guān)鍵技術(shù)。

1) 可靠的初至波識(shí)別與初至走時(shí)拾取方法技術(shù)。我們提出,弱初至掩埋在強(qiáng)噪聲的情形下,基于高維、多屬性、Markov最佳狀態(tài)轉(zhuǎn)移理論框架,充分利用機(jī)器學(xué)習(xí)領(lǐng)域發(fā)展出的新算法,實(shí)施初至波識(shí)別與走時(shí)拾取。不能認(rèn)為初至波FWI可以避開可靠的初至波識(shí)別與初至走時(shí)拾取的困難。事實(shí)上,復(fù)雜初至情形下初至波FWI也不能提供好的背景速度反演結(jié)果。除了初始模型差異過大之外,這正是初至波走時(shí)測量精度不能保證所引起的。單獨(dú)處理初至波識(shí)別與初至波走時(shí)拾取可以通過各種預(yù)處理和先驗(yàn)信息的引入解決實(shí)際數(shù)據(jù)(尤其是陸上地震數(shù)據(jù))不滿足初至波FWI理論假設(shè)引起的問題。到目前為止,弱初至掩埋在強(qiáng)噪聲情形下的初至波識(shí)別與初至走時(shí)拾取依然是沒有得到很好解決的問題。

2) 高效率且適應(yīng)非水平地表的初至波走時(shí)波動(dòng)理論層析成像方法技術(shù)。

3) (類)CMP道集動(dòng)校正掃描速度分析方法技術(shù)。類CMP道集是指PSTM或PSDM后成像道集反動(dòng)校正后得到的所謂的CMP道集,反動(dòng)校正用偏移速度。類CMP道集受橫向變速的影響小,動(dòng)校正速度分析結(jié)果更符合地質(zhì)實(shí)際。

圖2 基于特征波成像+寬帶波阻抗建模的定量地震波成像技術(shù)路線

4) 特征反射結(jié)構(gòu)導(dǎo)引下的特征反射波層析成像方法技術(shù)。利用反射波建立中深層背景速度不需要所有層位的反射波信息,基本原則是先用一些特征反射結(jié)構(gòu)(特征反射層)對應(yīng)的反射波同相軸(走時(shí)),進(jìn)行大套地層的層速度估計(jì)。然后,逐漸放入一些相對弱一些的反射層對應(yīng)的反射波同相軸(走時(shí)),提高層速度估計(jì)結(jié)果的分辨率。

5) 反射波FWI層析成像方法技術(shù)。在步驟4提供的初始速度模型的基礎(chǔ)上,進(jìn)行少許幾次反射波FWI層析成像迭代。其中的關(guān)鍵是要引入對梯度有地質(zhì)含義的正則化,并盡可能壓制掉偏移響應(yīng)梯度項(xiàng)的影響。

油氣地震勘探中,背景密度模型的構(gòu)建是比較困難的。原則上,選擇包含密度的波動(dòng)方程(譬如變密度聲波方程),利用FWI進(jìn)行多參數(shù)反演是可以反演密度的。但這是另一個(gè)反演騙局。實(shí)際數(shù)據(jù)情形下,密度變化會(huì)影響速度,速度影響波的到達(dá)時(shí),因而似乎可以由波的到達(dá)時(shí)反演密度。但是,密度和波的到達(dá)時(shí)的數(shù)學(xué)關(guān)系是什么?目前缺乏這樣的數(shù)學(xué)物理關(guān)系式。它們的關(guān)系可能是非線性的異常復(fù)雜的關(guān)系。密度變化影響同相軸上地震子波的振幅是確定的,密度變化導(dǎo)致波阻抗和反射子波振幅變化是確定的。但是,影響振幅的因素太多,這些因素構(gòu)成了由疊前地震數(shù)據(jù)基于Bayes估計(jì)理論下的FWI反演密度,因?yàn)榉蔷€性性很強(qiáng),極難在實(shí)用中見到效果的根本原因。據(jù)此,我們提出兩條構(gòu)建背景密度模型的技術(shù)路線。首先,由測井?dāng)?shù)據(jù)提供的密度信息和地震波成像提供的反射結(jié)構(gòu)信息,提出多信息約束下的密度場插值反問題,這是多約束條件下的非線性性優(yōu)化問題,解此優(yōu)化問題,得到背景密度場。其次,如果探區(qū)內(nèi)有重力勢場觀測數(shù)據(jù),再提一個(gè)重力勢場約束下多信息約束的密度反演問題,把前一步插值方法獲得的背景密度場作為初始值,求解該反問題得到所要的背景密度場。由背景速度場和背景密度場就可以得到背景波阻抗。為了方便,也可以借用已有的速度和密度之間的巖石物理關(guān)系式,計(jì)算出密度場,從而得到背景密度場。

寬帶反射系數(shù)(或高波數(shù)參數(shù)擾動(dòng)量)估計(jì)在寬帶波阻抗建模中提供高波數(shù)帶(譬如9~50Hz)的(相對)波阻抗變化信息。理論上,我們需要的是零角度反射系數(shù)。實(shí)際情況下,一般由角度反射系數(shù)的加權(quán)疊加作為零角度反射系數(shù)的近似。理論上,我們需要的是只有主瓣沒有旁瓣的零角度反射系數(shù),主瓣寬度也要盡可能地小。實(shí)際上,這也是地震波成像的期望目標(biāo)。有限頻帶的震源激發(fā)和有限孔徑的地表接收下的疊前地震數(shù)據(jù)不支持我們得到這樣的期望成像目標(biāo)。實(shí)際數(shù)據(jù)情況下,適用于寬帶波阻抗建模的反射系數(shù)獲取方法應(yīng)該是:①做好地表一致性處理;②做好保真的偏移成像(至少要消除球面擴(kuò)散和吸收衰減)得到最佳的Imig(保真的偏移成像結(jié)果);③用H*ITrue=Imig定義一個(gè)其它信息約束(主要是三維空間中的構(gòu)造約束)下的圖像反褶積問題(其中H是Hessian矩陣),求解該反問題得到更寬帶的Itrue(頻帶拓寬后的保真的偏移成像結(jié)果);④井震結(jié)合標(biāo)定Itrue的幅值大小。至此,可以認(rèn)為得到了定量的寬帶反射系數(shù)估計(jì)結(jié)果。

目前,油氣地震勘探地震地質(zhì)解釋的主要信息來源就是帶限反射系數(shù)成像結(jié)果。但是,很顯然反射系數(shù)R=(I2-I1)/(I2+I1)僅僅反映界面處的巖性變化。深度域全三維空間中大量的巖性變化信息并不包含在帶限反射系數(shù)成像結(jié)果中,而地質(zhì)解釋人員關(guān)心的是深度域全三維空間中介質(zhì)巖性的變化。這也是為什么一定要追求深度域全三維空間中巖性變化信息完整刻畫的地震波成像結(jié)果的根源。眾所周知,寬波數(shù)帶的縱波速度vP(x,y,z)、橫波速度vS(x,y,z)、密度ρ(x,y,z),波阻抗I(x,y,z),泊松比σ(x,y,z)都可以從不同的角度描述深度域全三維空間中的巖性變化,它們之間并不是相互獨(dú)立的,而是緊密關(guān)聯(lián)的。目前,我們的成像目標(biāo)是寬波數(shù)帶的波阻抗I(x,y,z)。事實(shí)上,寬波數(shù)帶的縱波速度vP(x,y,z)也是很合適的選擇。我們不是將它定為Bayes估計(jì)理論下FWI方法技術(shù)的成像目標(biāo),而是在已知背景阻抗和寬帶反射系數(shù)(或高波數(shù)參數(shù)擾動(dòng)量)估計(jì)結(jié)果的基礎(chǔ)上,提一個(gè)在其他先驗(yàn)信息約束下的信息融合反問題,進(jìn)一步融入更多的有地質(zhì)意義的先驗(yàn)信息(譬如地震相),將背景阻抗和寬帶反射系數(shù)融合成寬帶波阻抗。這樣做避開了Bayes估計(jì)理論下FWI方法技術(shù)估計(jì)寬帶波阻抗的強(qiáng)非線性性,得到的結(jié)果更穩(wěn)定、更符合地質(zhì)邏輯,而且沒有降低背景阻抗和寬帶反射系數(shù)所包含信息的分辨率。

3 寬帶波阻抗建模影響因素分析

定義并求解一個(gè)約束優(yōu)化問題實(shí)現(xiàn)信息融合寬帶波阻抗建模,目的是把背景阻抗和寬帶反射系數(shù)在其他信息的約束下,有機(jī)地融合在一起,得到低波數(shù)帶(0~3Hz)、中波數(shù)帶(4~8Hz)和高波數(shù)帶(9~50Hz)成分都合理存在的寬帶波阻抗模型。通過數(shù)值實(shí)驗(yàn)在Marmousi模型上分析背景阻抗的精度和寬帶反射系數(shù)的精度對寬帶波阻抗建模結(jié)果的影響。目的是依據(jù)分析結(jié)果得到提高寬帶波阻抗建模精度的合理做法。

3.1 背景阻抗的生成

將圖3所示的不同平滑半徑的高斯平滑算子施加在Marmousi模型上,得到不同平滑程度的背景速度。用描述速度密度之間關(guān)系的Gardner關(guān)系式產(chǎn)生密度場,從而生成平滑程度的背景阻抗。圖4展示了平滑半徑分別為40m,120m,240m,360m的背景阻抗。

圖3 半徑為R的Gauss平滑算子示意

為了更清楚地看出平滑算子的影響,從圖4中各子圖抽出第200道(模型中間一道)進(jìn)行阻抗曲線展示,其中圖5a為真值,圖5b至圖5e分別是40m,120m,240m和360m平滑后的結(jié)果。

由于頻率廣為大家熟知,通過關(guān)系式f/(df*nf)=k/(dk*nk)把深度域模型Fourier變換對應(yīng)的波數(shù)折算成頻率。其中,df,nf分別代表頻率域采樣間隔以及采樣點(diǎn)數(shù);dk,nk代表波數(shù)域采樣間隔以及采樣點(diǎn)數(shù)。圖6所示為圖5中各曲線對應(yīng)的振幅譜。

從圖6可以看出,隨著平滑尺度的增大,高波數(shù)的成分逐漸降低。平滑半徑R=40m時(shí),波阻抗中40Hz以上的頻率完全衰減;平滑半徑R=120m時(shí),16Hz以上的頻率完全衰減;平滑半徑R=240m時(shí)10Hz以上的頻率完全衰減;平滑半徑R=360m時(shí)4Hz以上的頻率基本完全衰減。

3.2 寬帶波阻抗建模需要什么樣的反射系數(shù)成像結(jié)果?

即便最高精度的疊前深度偏移成像,由于受實(shí)測地震數(shù)據(jù)中激發(fā)接收的地震子波頻帶有限和僅在地表觀測的觀測方式的孔徑有限,得到的帶限反射系數(shù)必然是雷克(Ricker)子波型反射系數(shù)。雷克(Ricker)子波型反射系數(shù)的典型特征是主瓣比較寬、旁瓣比較深。主要原因是缺低中頻和頻帶窄。我們認(rèn)為,各種地震子波的分辨率取決于完整的振幅譜,而相位譜的影響是次要的。低中頻優(yōu)勢成分的存在是讓子波旁瓣變淺,從而提高弱層分辨能力的重要信息源。由于雷克(Ricker)子波型反射系數(shù)主瓣比較寬、旁瓣比較深,因而它不是寬帶波阻抗建模所期望的反射系數(shù)。

圖4 不同平滑半徑的背景阻抗a 平滑尺度R=40m的阻抗結(jié)果; b 平滑尺度R=120m的阻抗結(jié)果; c 平滑尺度R=240m的阻抗結(jié)果; d 平滑尺度R=360m的阻抗結(jié)果

圖5 平滑后的背景阻抗抽線對比a 第200道的真值; b 40m平滑后的結(jié)果; c 120m平滑后的結(jié)果; d 240m平滑后的結(jié)果; e 360m平滑后的結(jié)果

圖6 圖5中各曲線對應(yīng)的振幅譜a 真實(shí)模型對應(yīng)的振幅譜(波數(shù)譜); b 平滑尺度R=40m的阻抗結(jié)果對應(yīng)的振幅譜(波數(shù)譜); c 平滑尺度R=120m的阻抗結(jié)果對應(yīng)的振幅譜(波數(shù)譜); d 平滑尺度R=240m的阻抗結(jié)果對應(yīng)的振幅譜(波數(shù)譜); e 平滑尺度R=360m的阻抗結(jié)果對應(yīng)的振幅譜(波數(shù)譜)

3.3 寬帶波阻抗建模需要什么樣的帶限反射系數(shù)?

為了回答這個(gè)問題,我們測試了4種類型的反射系數(shù)。這4種類型的反射系數(shù)見圖7。首先是脈沖型反射系數(shù),這是我們最想要得到的反射系數(shù),它是全頻帶的。其次是(我們WPI稱謂的)高斯型反射系數(shù),高斯型反射系數(shù)相對于脈沖型反射系數(shù)主要是衰減了高頻信息,但低中頻信息豐富。然后,是(我們WPI稱謂的)期望子波型反射系數(shù),具有窄主瓣、低旁瓣的特征。最后,是雷克(Ricker)子波型反射系數(shù)。這是寬帶波阻抗建模中最不想要的。雷克(Ricker)子波型反射系數(shù)的低頻和高頻都大幅度衰減,并且旁瓣的能量比較強(qiáng),會(huì)引入較多的假的波阻抗高頻信息。需要說明的是期望子波是我們根據(jù)地質(zhì)分辨率要求,用子波的主旁瓣寬度比和高度比為參數(shù)定義的高分辨子波。

圖7 四種類型的反射系數(shù)示意a 脈沖型反射系數(shù); b 高斯型反射系數(shù); c 期望子波型反射系數(shù); d 雷克子波型反射系數(shù)

圖8展示了上述4種反射系數(shù)表示的Marmousi模型。相當(dāng)于我們假定用疊前深度偏移成像方法得到了圖8展示的偏移成像剖面。由圖可見,脈沖型反射系數(shù)和高斯型反射系數(shù)表示的Marmousi模型成像剖面是理想的和期望的。

圖8 4種反射系數(shù)表示的Marmousi模型a 采用脈沖型反射系數(shù)表示的結(jié)果; b 采用高斯型反射系數(shù)表示的結(jié)果; c 采用期望子波型反射系數(shù)表示的結(jié)果; d 采用雷克子波型反射系數(shù)表示的結(jié)果

首先對4種反射系數(shù)剖面中CDP200道進(jìn)行抽線,并進(jìn)行頻譜分析,結(jié)果如圖9所示。可以明顯看出,脈沖型反射系數(shù)的頻帶展布較廣,全頻帶均有能量;高斯型反射系數(shù)在高頻端的能量明顯衰減,但低中頻成分有很好的保留;期望子波型反射系數(shù)在低頻和高頻端均有一定程度的衰減,但有效頻帶比較寬;雷克子波型反射系數(shù)的頻帶最窄,既缺低頻也缺高頻。

圖9 寬帶反射系數(shù)抽線及頻譜分析a 寬帶反射系數(shù)結(jié)果; b 寬帶反射系數(shù)頻譜(波數(shù)譜)

使用這4種反射系數(shù)剖面分別與平滑半徑為R=240m的背景阻抗進(jìn)行融合,得到的結(jié)果如圖10所示。此處,我們把信息融合寬帶波阻抗問題提成一個(gè)約束非線性優(yōu)化問題,具體提法如下:

(11a)

(11b)

圖10展示了不同寬帶反射系數(shù)與平滑半徑240m作用下的背景阻抗融合得到的阻抗結(jié)果;圖11展示了從CDP200處抽取的各種融合阻抗與對應(yīng)的背景阻抗和真阻抗的對比。從圖中可以看出,脈沖型反射系數(shù)與背景阻抗的融合結(jié)果趨于真解,與真波阻抗的相關(guān)度達(dá)到99.83%;高斯型反射系數(shù)由于高頻部分受損,融合阻抗的分辨率降低,但與真波阻抗的相關(guān)度達(dá)到99.47%;期望子波盡管是寬帶子波,但由于中低頻成分的缺失,旁瓣加深,融合結(jié)果中的低中頻成分難以補(bǔ)全,且已經(jīng)引入了假的高波數(shù)成分,與真波阻抗的相關(guān)度只能達(dá)到96.36%。雷克子波型反射系數(shù)融合效果最差,融合后波阻抗與真阻抗的相關(guān)度只有90.20%,從圖11中明顯可以看出它們之間的差異。

從融合過程可以清楚地看出,油氣地震勘探中寬帶成像子波對于介質(zhì)巖性刻畫的重要性,寬帶成像子波的獲取最重要的是既要有寬帶的激發(fā)接收子波,還要有寬方位長偏移距的觀測。尤其低中頻成分的缺失,成像結(jié)果對地下介質(zhì)巖性描述的誤差會(huì)很大。高密度和高信噪比的觀測對于巖性成像的貢獻(xiàn)是從其它角度體現(xiàn)的。

然后,再分析背景阻抗建模結(jié)果和寬帶反射系數(shù)成像結(jié)果對寬帶波阻抗建模結(jié)果的不同波數(shù)帶成分的貢獻(xiàn)。圖12展示了用不同類型寬帶反射系數(shù)剖面與背景阻抗(R=240m)融合后背景阻抗、融合阻抗和真阻抗的振幅譜對比。

圖10 不同寬帶反射系數(shù)與平滑半徑為240m的背景阻抗融合得到的阻抗結(jié)果a 采用脈沖型反射系數(shù)融合的波阻抗結(jié)果; b 采用高斯型反射系數(shù)融合的波阻抗結(jié)果; c 采用期望子波型反射系數(shù)融合的波阻抗結(jié)果; d 采用雷克子波型反射系數(shù)融合的波阻抗結(jié)果

圖11 不同寬帶反射系數(shù)融合結(jié)果抽線(CDP200)對比a 采用脈沖型反射系數(shù)融合的波阻抗結(jié)果抽線對比; b 采用高斯型反射系數(shù)融合的波阻抗結(jié)果抽線對比; c 采用期望子波型反射系數(shù)融合的波阻抗結(jié)果抽線對比; d 采用雷克子波型反射系數(shù)融合的波阻抗結(jié)果抽線對比

從圖12可以清楚地看出,融合阻抗的背景成分都是由背景阻抗貢獻(xiàn)出來的,其它部分是由反射系數(shù)成像結(jié)果貢獻(xiàn)出來的。這是很自然的邏輯結(jié)果。但是,中波數(shù)帶(譬如4~8Hz)也可以由反射系數(shù)成像補(bǔ)充出來,只要反射系數(shù)成像結(jié)果中包含有波阻抗中波數(shù)段的信息。事實(shí)上,脈沖型反射系數(shù)中包含了除零波數(shù)成分外所有波數(shù)成分的波阻抗變化信息。因此,道積分方法重建波阻抗時(shí),對于脈沖型反射系數(shù),理論上只要提供任意一層的正確波阻抗就可以得到真波阻抗的恢復(fù)。當(dāng)然,這又是一個(gè)反演騙局。實(shí)際上道積分方法受到的影響很多,根本不可能得到實(shí)際可用的波阻抗建模結(jié)果[14]。高斯型反射系數(shù)成像結(jié)果除了超高頻部分不能正確貢獻(xiàn)外,其他部分的波阻抗貢獻(xiàn)結(jié)果都是可靠的。應(yīng)該注意這個(gè)結(jié)論與高斯函數(shù)的方差是密切聯(lián)系的。較大方差下的高斯函數(shù)也會(huì)濾除過多的低中頻和超高頻成分,從而降低波阻抗建模的精度。理想子波型反射系數(shù)成像結(jié)果對寬帶阻抗的貢獻(xiàn),誤差較大的是中波數(shù)段。雷克子波型反射系數(shù)對寬帶阻抗的貢獻(xiàn),誤差較大的既包括中波數(shù)段,也包括高波數(shù)帶,高波數(shù)帶中明顯地引入了假的成分(圖11)。

上述融合過程再次告訴我們,地震波成像的真正目標(biāo)既包括精確的背景速度(背景阻抗),也包括最好是高斯型反射系數(shù)呈現(xiàn)出的成像結(jié)果(追求脈沖型反射系數(shù)成像結(jié)果過于理想)。即便得到高斯型反射系數(shù)呈現(xiàn)出的成像結(jié)果,成像分辨率也會(huì)有本質(zhì)性的提升。但是,高斯型反射系數(shù)成像結(jié)果也很難得到,疊前地震數(shù)據(jù)提供的信息不足以得到這樣的成像結(jié)果。無論如何,我們知道了地震波疊前偏移成像的目標(biāo)不是帶限的反射系數(shù),尤其不是具有寬主瓣、深(寬)旁瓣特征的雷克(Ricker)子波型的反射系數(shù),期望的偏移成像結(jié)果應(yīng)該是方差較小的高斯型反射系數(shù)。常規(guī)的偏移成像無論如何都不可能得到高斯型反射系數(shù)成像結(jié)果,但是約束反演成像方法有可能得到實(shí)用的高斯型反射系數(shù)成像結(jié)果,只不過要引入其它的先驗(yàn)信息。

3.3 背景阻抗對寬帶波阻抗建模的影響

為了分析不同平滑半徑獲得的背景阻抗(實(shí)質(zhì)上是不同精度的背景阻抗)對融合阻抗的精度影響,將圖4中4種背景阻抗與高斯型的寬帶反射系數(shù)(圖8)進(jìn)行融合,得到圖13所示的融合結(jié)果。圖14為CDP200處的抽線對比。

從圖13和圖14可以看出,4種不同平滑半徑得到的背景阻抗剖面與高斯型反射系數(shù)剖面的融合結(jié)果,精度都相當(dāng)高。說明當(dāng)?shù)卣鸪上衿拭姹憩F(xiàn)為高斯型反射系數(shù)時(shí),其中包含了波阻抗中的豐富的低中波數(shù)信息,即使用很平滑的背景阻抗,也能把寬波數(shù)帶的波阻抗重構(gòu)出來。各個(gè)波數(shù)帶中沒有顯著的成分缺失,重構(gòu)精度很高。

圖12 不同類型寬帶反射系數(shù)剖面與背景阻抗融合后,背景阻抗、融合阻抗和真阻抗的振幅譜對比a 對應(yīng)脈沖型反射系數(shù); b 對應(yīng)高斯型反射系數(shù); c 對應(yīng)理想子波型反射系數(shù); d 對應(yīng)雷克子波型反射系數(shù)

圖13 不同平滑尺度的背景阻抗與高斯型反射系數(shù)融合得到的融合阻抗a 采用平滑尺度R=40m的背景阻抗融合結(jié)果; b 采用平滑尺度R=120m的背景阻抗融合結(jié)果; c 采用平滑尺度R=240m的背景阻抗融合結(jié)果; d 采用平滑尺度R=360m的背景阻抗融合結(jié)果

圖14 CDP200處的抽線對比a 采用平滑尺度R=40m的背景阻抗融合結(jié)果抽線對比; b 采用平滑尺度R=120m的背景阻抗融合結(jié)果抽線對比; c 采用平滑尺度R=240m的背景阻抗融合結(jié)果抽線對比; d 采用平滑尺度R=360m的背景阻抗融合結(jié)果抽線對比

為了看到不同背景阻抗剖面與不同類型反射系數(shù)剖面融合結(jié)果之間的差異,又做了上述4種不同平滑半徑得到的背景阻抗剖面與期望子波型反射系數(shù)剖面(圖8)的融合,結(jié)果如圖15所示。圖16為CDP200處的抽線對比,圖17還對比了CDP200處各種融合阻抗、背景阻抗與真阻抗的振幅譜。

圖15給出的最主要的結(jié)論是:平滑半徑增大導(dǎo)致了背景波阻抗中中波數(shù)成分損失,理想子波型反射系數(shù)剖面也缺乏中波數(shù)成分,從而導(dǎo)致中波數(shù)阻抗不能被重構(gòu)(圖16和圖17)。中波數(shù)段成分的損失,使得重構(gòu)的波阻抗與真波阻抗出現(xiàn)明顯差異,重構(gòu)的波阻抗不能真實(shí)反映地下介質(zhì)巖性的變化。很顯然,拓寬成像子波的頻帶,提高背景阻抗的建模精度都是補(bǔ)充中波數(shù)阻抗成分的措施。

前述分析時(shí),背景阻抗成分是正確的。但是對于實(shí)際情況,背景阻抗也會(huì)存在一定的誤差。因此我們又分析了背景阻抗錯(cuò)誤情況下的信息融合寬帶波阻抗重構(gòu)結(jié)果。背景阻抗誤差引入方式為:真實(shí)阻抗進(jìn)行R=240m的高斯平滑后,整體減小15%;整體增加15%;加入15%量級(jí)的隨機(jī)擾動(dòng)。用高斯型反射系數(shù)成像結(jié)果與它們分別進(jìn)行信息融合重構(gòu),重構(gòu)結(jié)果如圖18和圖19所示。

當(dāng)背景阻抗為真實(shí)阻抗經(jīng)過半徑為240m的高斯平滑時(shí),重構(gòu)阻抗與真阻抗匹配度達(dá)到99.27%;但是若此時(shí)背景速度由上述背景速度整體減小15%時(shí),匹配度僅達(dá)到89.28%;若此時(shí)背景速度由上述背景速度整體增加15%時(shí),匹配度僅達(dá)到77.55%;若此時(shí)背景速度由上述背景速度加入15%量級(jí)的隨機(jī)擾動(dòng),對重構(gòu)結(jié)果影響不大,匹配度可達(dá)99.26%。值得慶幸的是,在實(shí)際情況下建立背景速度模型時(shí),出現(xiàn)系統(tǒng)誤差的概率不大。

圖15 不同平滑尺度的背景阻抗與期望子波型反射系數(shù)進(jìn)行融合得到的阻抗a 采用平滑尺度R=40m的背景阻抗融合結(jié)果; b 采用平滑尺度R=120m的背景阻抗融合結(jié)果; c 采用平滑尺度R=240m的背景阻抗融合結(jié)果; d 采用平滑尺度R=360m的背景阻抗融合結(jié)果

圖16 CDP200處的抽線對比a 采用平滑尺度R=40m的背景阻抗融合結(jié)果抽線對比; b 采用平滑尺度R=120m的背景阻抗融合結(jié)果抽線對比; c 采用平滑尺度R=240m的背景阻抗融合結(jié)果抽線對比; d 采用平滑尺度R=360m的背景阻抗融合結(jié)果抽線對比

圖17 CDP200處各種融合阻抗、背景阻抗與真阻抗的振幅譜a 平滑尺度R=40m時(shí)的對比情況; b 平滑尺度R=120m時(shí)的對比情況; c 平滑尺度R=240m時(shí)的對比情況; d 平滑尺度R=360m時(shí)的對比情況

圖18 錯(cuò)誤背景阻抗與高斯型反射系數(shù)進(jìn)行融合得到的阻抗a 采用真實(shí)波阻抗平滑240m后整體減小15%的背景波阻抗融合的結(jié)果; b 采用真實(shí)波阻抗平滑240m后的背景波阻抗融合的結(jié)果; c 采用真實(shí)波阻抗平滑240m后整體增大15%的背景波阻抗融合的結(jié)果; d 采用真實(shí)波阻抗平滑240m后整體加入15%量級(jí)的隨機(jī)擾動(dòng)的背景波阻抗融合的結(jié)果

圖19 錯(cuò)誤背景阻抗與高斯型反射系數(shù)進(jìn)行融合得到的阻抗在CDP200處的抽線對比a 采用真實(shí)波阻抗平滑240m后整體減小15%的背景波阻抗融合的結(jié)果抽線對比; b 采用真實(shí)波阻抗平滑240m后的背景波阻抗融合的結(jié)果抽線對比; c 采用真實(shí)波阻抗平滑240m后整體增大15%的背景波阻抗融合的結(jié)果抽線對比; d 采用真實(shí)波阻抗平滑240m后整體加入15%量級(jí)的隨機(jī)擾動(dòng)的背景波阻抗融合的結(jié)果抽線對比

4 寬帶波阻抗建模及定量地震波成像的基本認(rèn)識(shí)

油氣地震勘探的目的就是通過地震波成像提供的對地下介質(zhì)彈性參數(shù)的準(zhǔn)確刻畫,結(jié)合鉆井、測井、巖石物理和油氣地質(zhì)學(xué)知識(shí),準(zhǔn)確識(shí)別與評價(jià)油氣藏,進(jìn)行高效益的油氣開發(fā)。

深度域真三維空間中地下介質(zhì)彈性參數(shù)的準(zhǔn)確刻畫是地震波成像的本質(zhì)目標(biāo),也是油氣地震勘探中的核心問題。無論是深度域真三維空間中的寬波數(shù)帶縱波速度vP(x,y,z)或是寬帶波阻抗I(x,y,z)都能反映地下介質(zhì)巖性的變化,而帶限反射系數(shù)僅揭示介質(zhì)巖性在某些地方(波阻抗躍變的地方)發(fā)生突變,反映的是介質(zhì)波阻抗的相對變化量(不是絕對變化量)。很顯然,對于準(zhǔn)確識(shí)別與評價(jià)油氣藏而言,寬波數(shù)帶縱波速度vP(x,y,z)或是寬帶波阻抗I(x,y,z)建模結(jié)果比帶限反射系數(shù)更適合基礎(chǔ)數(shù)據(jù)。很長時(shí)間以來,地震波成像及后續(xù)的地震地質(zhì)解釋之所以總是基于帶限反射系數(shù)成像結(jié)果,首要原因是帶限反射系數(shù)的估計(jì)(或成像)是一個(gè)更穩(wěn)健的問題,更容易得到可靠的解。以帶限反射系數(shù)為目標(biāo)的偏移成像甚至不能稱為是一個(gè)真正的反問題。偏移成像是油氣地震勘探領(lǐng)域獨(dú)創(chuàng)的方法技術(shù)。遺憾的是,理論上說,它僅是一種地下繞(反)射點(diǎn)的空間定位方法技術(shù),振幅相對保真的成像結(jié)果只是偏移成像方法技術(shù)的附帶結(jié)果而已。

油氣地震勘探發(fā)展到現(xiàn)在,地震波成像進(jìn)入定量的、以估計(jì)介質(zhì)彈性參數(shù)絕對變化為目標(biāo)的階段已不可避免。Bayes估計(jì)理論下的FWI技術(shù)本質(zhì)上就是要達(dá)到此目標(biāo)。實(shí)際數(shù)據(jù)不滿足它的理論假設(shè),估計(jì)(全)寬波數(shù)帶的介質(zhì)彈性參數(shù)的絕對變化是一個(gè)強(qiáng)非線性問題,至少這兩個(gè)主要原因限制了它的實(shí)用化。我們提出的特征波成像(CWI)+寬帶波阻抗建模(WBIM)也是為了估計(jì)(全)寬波數(shù)帶的介質(zhì)彈性參數(shù)的絕對變化。當(dāng)前的目標(biāo)是信息融合寬帶波阻抗建模(或重構(gòu))。用石油工業(yè)的術(shù)語講,這是一種針對巖性儲(chǔ)層描述的、定量的地震波成像路線與方法技術(shù)系列。FWI的提法側(cè)重于強(qiáng)調(diào)算法;針對巖性儲(chǔ)層描述的、定量的地震波成像更強(qiáng)調(diào)油氣勘探的目標(biāo)和需求。

依據(jù)信息融合寬帶波阻抗建模(或重構(gòu))的思路,定量地震波成像被分成背景阻抗建模、寬帶反射系數(shù)估計(jì)和信息融合寬帶波阻抗建模3個(gè)核心環(huán)節(jié)。

第2節(jié)的分析指出,背景阻抗建模主要利用炮集中透射波(直達(dá)波、折射波和Diving Wave)和反(繞)射波的走時(shí)信息,求解線性反演問題(即所謂的層析成像),得到定量的背景速度估計(jì)?;谧邥r(shí)信息估計(jì)得到的背景速度與介質(zhì)層速度的背景部分在量級(jí)上是一致的,不是特別復(fù)雜的情形下,低波數(shù)的速度估計(jì)結(jié)果精度也是很高的。為了給背景速度補(bǔ)充進(jìn)更多的中波數(shù)成分,僅僅靠算法的改進(jìn)(譬如用FWI算法)和提高疊前數(shù)據(jù)的品質(zhì)(譬如“兩寬兩高”數(shù)據(jù))是不夠的,引入其他先驗(yàn)信息非常必要。目前認(rèn)為,最重要的先驗(yàn)信息就是地下介質(zhì)的結(jié)構(gòu)信息。首先是大套的反射層位信息、斷層信息和不規(guī)則地質(zhì)體信息,如果可能,再繼續(xù)引入地震相這樣更細(xì)節(jié)的地質(zhì)信息,會(huì)大幅度補(bǔ)充背景速度中的中波數(shù)段速度成分,從而提升整個(gè)背景速度建模的精度。上述分析可知,低中波數(shù)波阻抗成分的缺失對于定量地描述巖性變化會(huì)引起很大的誤差。在當(dāng)前“兩寬兩高”地震數(shù)據(jù)采集逐漸成為常態(tài),高精度的速度建模算法越來越多的情況下,通過機(jī)器學(xué)習(xí)算法提取更多的用于補(bǔ)充速度建模中波數(shù)成分的先驗(yàn)信息,進(jìn)而獲得可靠的、包含豐富的低、中波數(shù)信息的背景阻抗模型,可以保證最終定量的巖性成像結(jié)果的可靠性。

眾所周知,利用地震波同相軸上子波的振幅(波形)變化,估計(jì)介質(zhì)參數(shù)的擾動(dòng)是地震波成像中更關(guān)鍵的步驟。寬帶反射系數(shù)估計(jì)利用的就是炮集中反(繞)射波同相軸上地震子波的振幅信息。但是,影響同相軸上地震子波振幅的因素太多。理論上,應(yīng)該剝離其它因素對振幅的改變,僅留下介質(zhì)阻抗擾動(dòng)引起的振幅擾動(dòng),才能用適當(dāng)?shù)乃惴ㄓ刹▓稣穹鶖_動(dòng)反演地下介質(zhì)阻抗擾動(dòng)。

實(shí)際上,地震數(shù)據(jù)(波場)是經(jīng)過人工震源的激發(fā),大地介質(zhì)系統(tǒng)的濾波和檢波器系統(tǒng)的檢測和記錄才被獲取到的。這是一個(gè)非常復(fù)雜的過程。通過對震源和檢波器系統(tǒng)的簡化,把地震波成像(或參數(shù)反演)更廣義地抽象如圖20所示。

觀測的地震波場可以抽象表示為:

d(xR,t,xS)=C[u(x,t,xS)]+η

(12)

A(m)u=f

(13)

式中:f是震源函數(shù),u為波動(dòng)方程的解。A為人為構(gòu)建的波動(dòng)方程,C為地震數(shù)據(jù)的采樣算子,將連續(xù)的地震波場轉(zhuǎn)換為指定觀測點(diǎn)的離散的空間時(shí)間函數(shù),η為采集過程引入的噪聲函數(shù)。

圖20 反演成像系統(tǒng)示意

忽略對震源激發(fā)的振幅的影響分析,通過對檢波器系統(tǒng)的抽象過程可以看出:

1) 對地面振動(dòng)位移的感知,檢波器將地面震動(dòng)位移線性地轉(zhuǎn)化為電信號(hào)。在該過程中,要求電信號(hào)強(qiáng)度與振動(dòng)強(qiáng)度存在常比例關(guān)系,并且需要檢波器響應(yīng)動(dòng)態(tài)范圍較大,最強(qiáng)振動(dòng)與最弱振動(dòng)都能夠被有效地感知到并按上述要求產(chǎn)生電信號(hào)。此時(shí)需要注意的是,地震子波(波前面)平穩(wěn)下來之后才能把地面的近似振動(dòng)幅度當(dāng)成源的強(qiáng)度,平穩(wěn)之前的振動(dòng)位移應(yīng)該舍棄。因此,能否保真感知并產(chǎn)生電信號(hào),動(dòng)態(tài)范圍是否能適應(yīng)勘探需求,弱振動(dòng)能否保真地記錄下來,都使檢波器的感知面臨巨大的挑戰(zhàn)。

2) 能否把檢波器輸出的被認(rèn)為保真感知的電信號(hào)轉(zhuǎn)換并記錄為有效數(shù)字,最關(guān)鍵的是模擬轉(zhuǎn)換設(shè)備。在當(dāng)前高性能計(jì)算機(jī)系統(tǒng)下,模擬轉(zhuǎn)換數(shù)字信號(hào)的存儲(chǔ)并不存在問題。然而模擬轉(zhuǎn)換的動(dòng)態(tài)范圍是個(gè)問題,過粗的采樣忽略了弱信號(hào)。也說明感知到的大動(dòng)態(tài)范圍的振動(dòng)(電信號(hào))不一定能被保真地轉(zhuǎn)化成數(shù)字信號(hào)記錄下來。因此,在檢波器的設(shè)計(jì)過程中,需要考慮深層弱信號(hào)的動(dòng)態(tài)范圍、感知部分能做到保真感知的動(dòng)態(tài)范圍和記錄部分能做到保真記錄的動(dòng)態(tài)范圍。

另外,檢波器的方向選擇性(不同分量對振動(dòng)位移方向的感知不同)能否應(yīng)對野外觀測的實(shí)際變化和檢波器是否能與地表介質(zhì)耦合以及引起的振幅畸變會(huì)有多大,對保真成像帶來了巨大的挑戰(zhàn)。

定量的地震波成像不仔細(xì)考察震源激發(fā)與檢波器接收的振幅保真性是不可能做好的。

3) 基于線性成像系統(tǒng)分析保真成像的影響因素。首先,分析觀測系統(tǒng)中的單個(gè)炮檢對情況下成像系統(tǒng)的照明及分辨率。在高頻近似下,單個(gè)炮檢對對應(yīng)的射線與地下某一點(diǎn)作用的幾何關(guān)系如圖21所示。

圖21 高頻近似下,單個(gè)炮檢對對應(yīng)的射線與地下單點(diǎn)作用的幾何關(guān)系

(14)

式中:ISR代表照明矢量,nSR為照明矢量方向上的單位矢量,v0(x)為散射點(diǎn)處的背景速度。照明矢量模的長短決定著分辨率的高低。可以看出,散射張角的大小對照明矢量的模長有重要影響,對于相同的背景速度,當(dāng)Θ=0時(shí),照明矢量的模最長,分辨率最高,此時(shí)對應(yīng)正入射(自激自收)情形;當(dāng)散射張角增大,照明矢量的模變短,分辨率降低;當(dāng)Θ=180°時(shí),照明矢量長度為零,此時(shí)對應(yīng)透射情形,對散射點(diǎn)沒有成像分辨率。散射張角對分辨率的影響可以解釋角度域成像道集(或偏移距域成像道集)中,大角度(或大偏移距)產(chǎn)生子波拉伸。另外,背景速度的大小也影響著照明矢量模的長度,背景速度越大,照明矢量模長越短,分辨率越低。一般情況下,地下介質(zhì)從淺到深,速度逐漸增大,解釋了通常深層成像結(jié)果分辨率比淺層低的原因。

觀測系統(tǒng)中,通過一對炮檢對在地下某一點(diǎn)對應(yīng)的照明矢量ISR,可以計(jì)算出該炮檢對在地下該點(diǎn)對應(yīng)的波數(shù)矢量kSR,即

kSR=ωISR

(15)

式中:ω代表角頻率。圖22為單個(gè)炮檢對波數(shù)矢量的投影。通過上式,可將地震波場波前面上(時(shí)間域)的子波的頻譜投影至波數(shù)域,波數(shù)矢量的方向與照明矢量相同,其幅值取決于波場子波的頻譜以及照明矢量的模長。波數(shù)矢量是分辨率分析中的關(guān)鍵量,將波數(shù)域中的kSR進(jìn)行Fourier反變換,可在空間域成像點(diǎn)處,得到一個(gè)帶限局部平面波前面。波前面延續(xù)方向與雙程走時(shí)場T(xS,xR,x)的等時(shí)面(偏移畫弧曲面)在成像點(diǎn)位置相切,與波數(shù)矢量kSR的方向正交。

圖22 單個(gè)炮檢對波數(shù)矢量的投影a 成像點(diǎn)處,地震波場波前面上時(shí)間域地震子波的頻譜(此處為主頻30Hz Ricker子波的頻譜); b 假定成像點(diǎn)處速度為2000m/s,炮檢對對應(yīng)的散射張角為0,利用照明矢量以及投影關(guān)系,將a中的頻譜投影至波數(shù)域,獲得的波數(shù)矢量,注意波數(shù)矢量展布的方向與選定的照明矢量方向相同; c 通過Fourier反變換,將b中波數(shù)矢量變換到空間域,得到帶限的局部平面波波前面,波前面的“厚度”決定了成像分辨率,其與照明矢量模的長度以及子波頻帶有關(guān),波前面的延續(xù)方向與波數(shù)矢量方向正交

由上述分析可知,子波頻帶、成像點(diǎn)處速度以及散射張角是影響分辨率的3個(gè)重要因素。

實(shí)際情形下,由于觀測系統(tǒng)的孔徑有限,并且炮檢點(diǎn)分布不均勻,加上地下復(fù)雜的速度結(jié)構(gòu)對波傳播路徑的影響,理想照明情形無法實(shí)現(xiàn)。實(shí)際上,觀測系統(tǒng)中在地下每個(gè)成像點(diǎn)處滿足某一反射角及方位角的炮檢對,對成像點(diǎn)的角度照明是不完備且不均勻的。圖23展示了在僅有一定角度照明范圍的情形下,波數(shù)矢量的展布及其對應(yīng)的空間域點(diǎn)擴(kuò)散函數(shù)(point spread function,PSF)。從圖中可以看出,將所有波數(shù)矢量進(jìn)行矢量求和,可以得到一個(gè)平均波數(shù)矢量kmean,其方向決定了空間PSF的高分辨率軸以及低分辨率軸。高分辨率軸沿著kmean的方向,在該方向上,PSF具有最高的分辨率,PSF在該方向上的“厚度”決定了成像系統(tǒng)的“縱向”分辨率,“縱向”分辨率受到子波頻帶以及照明矢量模長的影響;低分辨率軸的方向垂直于kmean的方向,在該方向上PSF的分辨率最低,PSF在該方向上的延續(xù)寬度決定了成像系統(tǒng)的“橫向”分辨率,“橫向”分辨率受到波數(shù)矢量在角度上覆蓋范圍的影響,角度覆蓋范圍越大,“橫向”分辨率越高。由圖23可以看出,所謂的“縱向”與“橫向”,實(shí)質(zhì)上是沿kmean以及垂直kmean的方向。

從上述分析可知,影響成像系統(tǒng)在某一成像點(diǎn)分辨率的本質(zhì)原因,是觀測系統(tǒng)中的炮檢對在該成像點(diǎn)對應(yīng)的波數(shù)矢量在波數(shù)域的展布范圍。展布范圍包含兩方面,一是單個(gè)波數(shù)矢量的長度,其次是所有波數(shù)矢量在不同角度方向上的分布范圍。成像點(diǎn)處的速度越小、散射張角越小、成像點(diǎn)處波前面上子波的頻帶越寬時(shí),波數(shù)矢量的長度越長,對應(yīng)的PSF有越高的“縱向”分辨率。其中,波前面上的子波頻帶取決于震源子波頻帶以及波傳播過程中的吸收衰減、頻散等效應(yīng)。觀測系統(tǒng)中的炮檢對,在成像點(diǎn)處對應(yīng)的波數(shù)矢量在波數(shù)域分布的角度范圍越寬,對應(yīng)的PSF有越高的“橫向”分辨率。波數(shù)矢量的展布方向取決于炮檢對相對于地下成像點(diǎn)的位置,同時(shí)也受到背景速度對入射、散射波場傳播方向的影響。綜上所述,從觀測系統(tǒng)設(shè)計(jì)的角度而言,想要使成像系統(tǒng)有較高的分辨率以及較均衡的照明,震源子波需要有寬頻帶,炮檢對的分布需要有較寬的孔徑,并且排列要盡可能規(guī)則。另外,為了使得成像結(jié)果有較高的信噪比,需要高密度的炮檢點(diǎn)的分布。對于高分辨率的成像系統(tǒng),需要“兩寬一高”(寬頻帶、寬方位、高密度)的觀測系統(tǒng)。

圖23 在不完備的角度照明情形下,波數(shù)矢量在波數(shù)域的展布及其對應(yīng)的空間PSF

綜上所述,由于震源函數(shù)未知(實(shí)際采集過程中震源函數(shù)變化較為復(fù)雜)、聲波波動(dòng)方程過于簡單(地震波在復(fù)雜介質(zhì)中的傳播需要復(fù)雜的波動(dòng)方程描述),檢波器系統(tǒng)的振幅保真問題,導(dǎo)致觀測數(shù)據(jù)中的數(shù)據(jù)誤差η遠(yuǎn)偏離高斯分布,使得保真成像變得異常困難。

我們認(rèn)為,寬帶反射系數(shù)估計(jì)不是一個(gè)數(shù)學(xué)算法可以解決的,需要一個(gè)合理的技術(shù)流程,結(jié)合合理的關(guān)鍵技術(shù),基本上能得到與地下介質(zhì)波阻抗擾動(dòng)量一致的、定量的成像結(jié)果。為此提出如下的寬帶反射系數(shù)定量化估計(jì)的技術(shù)流程。

1) 地表一致性振幅校正。消除激發(fā)端和接收端導(dǎo)致的地震子波的振幅變化,余下的振幅變化僅由傳播過程和阻抗擾動(dòng)引起;

2) 球面擴(kuò)散校正。偏移成像過程能處理球面擴(kuò)散校正時(shí),不必單獨(dú)處理;

3) 吸收衰減校正。偏移成像過程能處理吸收衰減校正時(shí),不必單獨(dú)處理;

4) 透射損失校正、薄層干涉校正和散射損失校正。這些因素引起的振幅變化往往被忽略;

5) 保真的疊前深度偏移成像。當(dāng)前,我們的建議是做好常規(guī)的疊前深度偏移(PSDM),產(chǎn)生(方位)角度道集,進(jìn)行照明優(yōu)選的保真疊加。不建議用LS_PSDM和LS_RTM,因?yàn)檫@二種方法施加的理論假設(shè)過多,實(shí)際數(shù)據(jù)很難滿足;

6) 成像域點(diǎn)擴(kuò)散函數(shù)(PSF)反褶積。目的是校正照明不均勻,進(jìn)一步提升反射系數(shù)成像結(jié)果保真度和分辨率;

7) 進(jìn)行測井反射系數(shù)與地震成像反射系數(shù)的量級(jí)標(biāo)定。我們建議在上述保真成像盡可能做好的情況下,僅僅對大套標(biāo)志性反射層的反射系數(shù)進(jìn)行一致性標(biāo)定。

總之,帶限反射系數(shù)的相對保真性由保真成像處理各環(huán)節(jié)保障。帶限反射系數(shù)的頻帶展寬,甚至擴(kuò)展成高斯型反射系數(shù),需要提一個(gè)合適的約束非線性反問題并進(jìn)行有效的求解;寬帶反射系數(shù)成像結(jié)果的定量化要靠測井反射系數(shù)標(biāo)定。可以認(rèn)為,井反射系數(shù)是一種質(zhì)控信息。

在上述認(rèn)識(shí)基礎(chǔ)上,我們將特征波成像(CWI)+寬帶波阻抗建模(WBIM)流程及技術(shù)系列應(yīng)用于理論和實(shí)際數(shù)據(jù)的處理,驗(yàn)證了方法技術(shù)系列的有效性。

5 定量地震波成像數(shù)值實(shí)驗(yàn)

我們對西北某探區(qū)抽象出的三維模型數(shù)據(jù)和三維實(shí)際數(shù)據(jù),用特征波成像(CWI)+寬帶波阻抗建模(WBIM)技術(shù)流程進(jìn)行數(shù)據(jù)處理。

5.1 三維模型數(shù)據(jù)寬帶波阻抗建模

首先,利用三維模型來模擬西北某探區(qū)實(shí)際地質(zhì)情況。其中,背景速度由模型真速度經(jīng)過半徑為375m的高斯平滑產(chǎn)生,用Garder速度密度關(guān)系生成背景密度,從而得到背景阻抗。寬帶反射系數(shù)rmig由Kirchhoff疊前深度偏移結(jié)果生成。利用(11)式定義的信息融合寬帶波阻抗重構(gòu)方法生成寬帶波阻抗建模結(jié)果。

圖24顯示了三維模型數(shù)據(jù)寬帶波阻抗建模的結(jié)果。可以看出,重構(gòu)阻抗和真實(shí)阻抗有良好的一致性,圖25的抽線對比了寬帶波阻抗建模的結(jié)果。由圖25可以看出,重構(gòu)的寬帶阻抗在低波數(shù)段和高波數(shù)段都是比較可靠的,5~10Hz左右的中波數(shù)段重構(gòu)阻抗的振幅譜值與真實(shí)阻抗的振幅譜值差異比較大??梢哉J(rèn)定,這是因?yàn)镵irchhoff疊前深度偏移成像結(jié)果中沒有提供中波段的信息。其中,由于帶限反射系數(shù)成像結(jié)果頻帶有限,反射系數(shù)旁瓣太深,導(dǎo)致波阻抗建模中出現(xiàn)假的波阻抗躍變。在之前的影響因素分析中已經(jīng)介紹過。中波數(shù)段部分缺失也可以歸結(jié)為半徑為375m的平滑半徑過大,導(dǎo)致背景阻抗中也沒有包含這部分中波數(shù)成分的信息。

圖24 三維模型數(shù)據(jù)寬帶波阻抗建模結(jié)果的剖面顯示a 反射系數(shù)振幅標(biāo)定后的Kirchhoff疊前深度偏移結(jié)果剖面; b 基于真速度平滑(半徑375m)后的背景速度產(chǎn)生的背景阻抗; c 信息融合寬帶阻抗建模結(jié)果; d 真阻抗模型

由圖26可以看出,振幅標(biāo)定后的Kirchhoff疊前深度偏移成像結(jié)果的深度切片并沒有很好地展示出地下巖性變化的真實(shí)情況。與真實(shí)阻抗相比,即便在幾何形態(tài)上,也與真實(shí)的阻抗擾動(dòng)差異較大。尤其在7.35km的深度切片上,這一點(diǎn)看得很清楚。理論上,這是不應(yīng)該出現(xiàn)的現(xiàn)象。如果疊前深度偏移做得符合理論預(yù)期,疊前深度偏移結(jié)果至少在幾何形態(tài)上能很準(zhǔn)確地反映地下介質(zhì)波阻抗擾動(dòng)情況。問題當(dāng)然出在Kirchhoff疊前深度偏移成像受偏移速度和成像算法的影響上。應(yīng)該提及的是,由于三維模型正演模擬數(shù)據(jù)量巨大,此處用的Kirchhoff疊前深度偏移成像結(jié)果是由模擬數(shù)據(jù)一方提供的。即便如此,寬帶波阻抗重構(gòu)結(jié)果與真實(shí)阻抗的一致性很高,顯然更適于巖性油藏的地質(zhì)解釋。

圖25 對應(yīng)圖24d虛線所示位置處的寬帶波阻抗建模結(jié)果的抽線對比a 圖24d左邊虛線位置處背景阻抗、重構(gòu)阻抗和真阻抗的對比; b 背景阻抗、重構(gòu)阻抗和真阻抗的振幅譜對比; c 圖24d中間虛線位置處背景阻抗、重構(gòu)阻抗和真阻抗的對比; d 背景阻抗、重構(gòu)阻抗和真阻抗的振幅譜對比; e 圖24d右邊虛線位置處背景阻抗、重構(gòu)阻抗和真阻抗的對比; f 背景阻抗、重構(gòu)阻抗和真阻抗的振幅譜對比

圖26 三維模型數(shù)據(jù)寬帶波阻抗建模結(jié)果不同深度的切片顯示a 7.35km深處標(biāo)定后的偏移成像結(jié)果顯示; b 7.35km深處重構(gòu)寬帶波阻抗切片顯示; c 7.35km深處真波阻抗切片顯示; d 8.675km深處標(biāo)定后的偏移成像結(jié)果顯示; e 8.675km深處重構(gòu)寬帶波阻抗切片顯示; f 8.675km深處真波阻抗切片顯示

這就從理論模型數(shù)據(jù)上證明了我們提出的信息融合寬帶波阻抗建模路線及方法技術(shù)的有效性。

5.2 實(shí)際數(shù)據(jù)寬帶波阻抗建模

選擇西北某探區(qū)的實(shí)際數(shù)據(jù)進(jìn)行寬帶波阻抗建模處理。關(guān)鍵技術(shù)環(huán)節(jié)包括:①背景波阻抗模型的建立;②波動(dòng)方程疊前深度偏移;③方位角度成像道集的輸出;④成像道集后處理,包括道集拉平、最佳照明選擇、最優(yōu)加權(quán)疊加;⑤成像域反褶積;⑥基于信息融合的寬帶波阻抗建模。

背景阻抗建模中的背景速度來自于疊前深度偏移所用的速度場,背景密度用Gardner關(guān)系由背景速度轉(zhuǎn)換得到。寬帶反射系數(shù)由前段中的步驟②~⑤產(chǎn)生。步驟⑥由求解(11)式完成。

圖27按剖面展示了實(shí)際數(shù)據(jù)信息融合寬帶波阻抗重構(gòu)結(jié)果。圖28以深度切片的方式展示了實(shí)際數(shù)據(jù)信息融合寬帶波阻抗重構(gòu)結(jié)果。由于不知道真實(shí)阻抗,信息融合寬帶波阻抗重構(gòu)結(jié)果好壞的判斷比較困難,只能從成像分辨率和地質(zhì)邏輯合理性來評判。從信息融合的角度,算法不會(huì)傷害原有信息的分辨率。理論上,反射系數(shù)是絕對阻抗的微分,反射系數(shù)剖面的視覺分辨率應(yīng)該會(huì)高于寬帶阻抗。但是,從圖28 展示的深度切片看,反射系數(shù)成像結(jié)果的地質(zhì)合理性不如重構(gòu)阻抗,在5.5km處顯得十分雜亂。寬帶阻抗成像結(jié)果的地質(zhì)合理性還要與該探區(qū)的地質(zhì)解釋人員結(jié)合再做進(jìn)一步判斷。

從實(shí)際數(shù)據(jù)處理結(jié)果看,我們提出的信息融合寬帶波阻抗建模路線及方法技術(shù)是可行的。但寬帶阻抗成像結(jié)果在油藏描述中的有效性,還要進(jìn)一步檢驗(yàn)。

無論如何,實(shí)際數(shù)據(jù)情況下,針對巖性儲(chǔ)層的定量地震波成像,要求地震數(shù)據(jù)采集環(huán)節(jié)、地震數(shù)據(jù)預(yù)處理各環(huán)節(jié)、背景速度建模、背景密度建模、寬帶反射系數(shù)定量估計(jì),都要全面提升質(zhì)量和精度,才有可能得到能顯著提升油藏描述精度和勘探效益的結(jié)果。因此,沿著針對巖性儲(chǔ)層的定量地震波成像方向提升勘探技術(shù),比用FWI技術(shù)引領(lǐng)勘探技術(shù)的發(fā)展更合理。

圖27 實(shí)際數(shù)據(jù)信息融合寬帶波阻抗重構(gòu)結(jié)果a 背景阻抗剖面; b 標(biāo)定后的寬帶反射系數(shù)剖面; c 重構(gòu)的寬帶波阻抗剖面; d 背景阻抗剖面; e 標(biāo)定后的寬帶反射系數(shù)剖面; f 重構(gòu)的寬帶波阻抗剖面

圖28 實(shí)際數(shù)據(jù)信息融合寬帶波阻抗重構(gòu)結(jié)果的深度切片顯示a 振幅標(biāo)定后的疊前深度偏移結(jié)果深度切片(5.5km); b 重構(gòu)波阻抗深度切片(5.5km); c 振幅標(biāo)定后的疊前深度偏移結(jié)果深度切片(7.9km); d 重構(gòu)波阻抗深度切片(7.9km); e 振幅標(biāo)定后的疊前深度偏移結(jié)果深度切片(8.1km); f 重構(gòu)波阻抗深度切片(8.1km)

6 結(jié)論與討論

隨著油氣勘探的新領(lǐng)域(深水、深層、復(fù)雜儲(chǔ)層、非常規(guī)儲(chǔ)層等)逐步擴(kuò)展,油氣勘探的目標(biāo)越來越復(fù)雜。但是,精確描述油氣藏、確定最佳井位和最佳鉆井方案、獲得最高油氣勘探效益的需求并沒有變化;精確描述油氣藏需要更精確的地震波成像結(jié)果的需求也沒有改變??傮w上講,提供寬波數(shù)帶的、能反映介質(zhì)巖性變化的、定量的地震波成像結(jié)果是油氣勘探的核心問題。

地震波成像的本質(zhì)就是利用疊前地震數(shù)據(jù)和與待估計(jì)參數(shù)相關(guān)的先驗(yàn)信息,利用人為提出的地下介質(zhì)參數(shù)變化和實(shí)測數(shù)據(jù)中波場變化之間的數(shù)學(xué)物理關(guān)系,利用Bayes參數(shù)估計(jì)理論下提出的(非)線性化反演方法,得到介質(zhì)參數(shù)變化的估計(jì)結(jié)果。FWI方法、層析成像方法、LS_PSDM/LS_RTM和一般的偏移成像就是典型的地震波成像方法?;跁r(shí)距關(guān)系的動(dòng)校正和靜校正本質(zhì)上也應(yīng)歸于地震波成像的范疇。寬波數(shù)帶的參數(shù)估計(jì)結(jié)果是地震波成像的根本目標(biāo)。帶限反射系數(shù)成像在石油工業(yè)界廣為應(yīng)用,主要是因?yàn)閷?yīng)的偏移成像算法穩(wěn)定而且物理意義清楚,并非帶限反射系數(shù)成像結(jié)果更能滿足精確油藏描述的需求,僅是深度域真三維空間中波阻抗相對變化的體現(xiàn)。實(shí)際上,在一個(gè)三維油氣勘探空間中,僅有很少部位有顯著的波阻抗相對變化。當(dāng)然,油氣勘探主要是對波阻抗發(fā)生變化的部位感興趣。

寬帶的、定量的地震波成像結(jié)果(譬如縱波速度vP(x,y,z)或波阻抗I(x,y,z))是在深度域真三維空間中對地下介質(zhì)巖性變化進(jìn)行準(zhǔn)確刻畫的物理量,具有明確的油氣地質(zhì)意義,適于進(jìn)行高精度的、定量的油藏描述。FWI反演方法的真正目的也就是獲得寬帶的、定量的參數(shù)估計(jì)結(jié)果,但這是個(gè)強(qiáng)非線性反問題,極難得到穩(wěn)定收斂的、有地質(zhì)意義的、能實(shí)用的解。

李慶忠[14]在《走向精確勘探的道路》一書中詳細(xì)分析了波阻抗建模的各種實(shí)際影響因素,很受啟發(fā)。但他的工作還沒有提升到信息融合的高度,沒有奠基在FWI反演成像理論基礎(chǔ)上,留下了很多未解的疑問。

我們提出的特征波成像(CWI)+信息融合寬帶波阻抗建模(WBIM)方法技術(shù)系列,本質(zhì)上是把FWI寬帶波阻抗反演,分解成背景速度建模、背景密度建模、線性化的寬帶反射系數(shù)估計(jì)這3個(gè)凸性更好的反問題進(jìn)行求解??紤]到實(shí)際數(shù)據(jù)的復(fù)雜性,這3個(gè)關(guān)鍵技術(shù)中的任何一個(gè)技術(shù)都還需要做進(jìn)一步的變通,目的是在實(shí)際數(shù)據(jù)的限制下還能得到“最佳”的成像結(jié)果。關(guān)鍵是我們并沒有把寬帶波阻抗估計(jì)提成一個(gè)非線性的反演問題來求解,而是轉(zhuǎn)化成一個(gè)信息融合問題。將已有的、認(rèn)為最好的背景速度建模、背景密度建模、線性化的寬帶反射系數(shù)估計(jì)結(jié)果,通過提出和求解一個(gè)約束非線性優(yōu)化問題,重構(gòu)出寬帶波阻抗模型。這樣做完全避開了反演求解寬帶波阻抗時(shí)的不穩(wěn)定和不收斂問題,同時(shí)在此過程中還可以進(jìn)一步引入其他先驗(yàn)信息,在信息融合這一步繼續(xù)提升寬帶阻抗建模的精度。

本文提出的由背景速度、背景密度和寬帶反射系數(shù),結(jié)合其他先驗(yàn)信息,融合后得到寬帶波阻抗,是不同于目前常規(guī)波阻抗反演的一種新的寬帶波阻抗建模方法技術(shù)。建模結(jié)果穩(wěn)定,地質(zhì)含義更清楚。寬帶波阻抗模型各波數(shù)帶的信息來源清楚,低頻(低波數(shù))來自背景波阻抗,高頻(高波數(shù))來自高分辨反射系數(shù)成像。中頻成分的補(bǔ)充依賴于背景速度建模和構(gòu)造成像精度的提高(背景速度的波數(shù)帶向中波數(shù)段展寬),也依賴于高分辨反射系數(shù)成像精度的提升(寬帶反射系數(shù)的波數(shù)帶向中波數(shù)段展寬)。

期望的寬帶反射系數(shù)是高斯型反射系數(shù),其次是理想子波型寬帶反射系數(shù)。Ricker子波型反射系數(shù)成像結(jié)果因?yàn)橛休^寬的主瓣和較為寬深的旁瓣,導(dǎo)致在寬帶波阻抗建模中引入很大的誤差。Ricker子波型反射系數(shù)成像結(jié)果之所以有這樣的特點(diǎn)主要是當(dāng)前很多實(shí)際觀測數(shù)據(jù)缺低中頻(譬如2~8Hz)有效成分,且觀測方位不夠?qū)捄推凭嗖粔蜷L導(dǎo)致的。數(shù)據(jù)的低信噪比對成像精度的影響首先體現(xiàn)在影響速度建模的精度上,然后體現(xiàn)在影響偏移成像的精度上。高精度寬帶波阻抗建模需要有低中頻(中波數(shù))段優(yōu)勢成分的寬帶反射系數(shù)(譬如方差較小的高斯型反射系數(shù))。寬帶反射系數(shù)中的低頻(包括中頻)成分的來源與寬頻子波、寬方位長偏移觀測又是密切關(guān)聯(lián)的。背景速度的估計(jì)精度事實(shí)上也嚴(yán)重依賴寬頻子波和寬方位長偏移觀測?!皟蓪拑筛摺钡卣饠?shù)據(jù)觀測對定量地震波成像具有核心作用。另外,即便做了保真的地震波成像,寬帶反射系數(shù)(量級(jí))依然需要測井反射系數(shù)的標(biāo)定,寬帶反射系數(shù)的量級(jí)對波阻抗反演結(jié)果影響很大。說明定量的地震波反演成像需要井與震的結(jié)合,這是無法回避的。

盡管“兩寬兩高”地震數(shù)據(jù)采集已經(jīng)成為業(yè)界的共識(shí),也已經(jīng)得到廣泛應(yīng)用,但是,滿足油藏描述期望的、定量的寬帶成像結(jié)果依然很難達(dá)成。我們一直認(rèn)為,寬帶地震波成像(或?qū)拵?shù)估計(jì))是(一直是)一個(gè)信息不足情況下的最佳估計(jì)問題。高精度的地震波成像(寬帶參數(shù)估計(jì))必須要引入更多的地質(zhì)、測井和巖石物理信息。針對巖性儲(chǔ)層描述的、定量的地震波成像過程與油藏描述過程逐漸融為一體。在大數(shù)據(jù)、人工智能和機(jī)器學(xué)習(xí)算法的賦能下,測井?dāng)?shù)據(jù)中與巖性變化相關(guān)的信息要充分地挖掘出來;地震成像剖面上與地下構(gòu)造、沉積相關(guān)的信息要盡可能地挖掘出來;巖石物理知識(shí),也包括定性的地質(zhì)知識(shí)也要能被計(jì)算機(jī)表達(dá),所有這些信息都會(huì)有助于提升寬帶參數(shù)估計(jì)的精度,在深度域全三維空間中盡可能準(zhǔn)確地描述地下介質(zhì)巖性的狀態(tài)。因此,今后的地震波成像技術(shù)發(fā)展,我們認(rèn)為應(yīng)該在定量的巖性成像概念引領(lǐng)下向前發(fā)展,而不是在FWI的理念下向前發(fā)展。追求定量的巖性成像,針對的是石油工業(yè)的需求;追求FWI成像,更多的是算法方面的深究。

盡管文中數(shù)值實(shí)驗(yàn)顯示我們提出的特征波成像(CWI)+信息融合寬帶波阻抗建模(BWIM)這種定量的巖性成像方法技術(shù)路線是合理且有效的,但要實(shí)現(xiàn)在油氣勘探中廣泛應(yīng)用的目標(biāo),還有很多工作要做。我們認(rèn)為當(dāng)前要深化的工作主要不是體現(xiàn)在發(fā)展更復(fù)雜的反演成像算法上,而是充分利用機(jī)器學(xué)習(xí)算法提取更多的能提升寬帶參數(shù)估計(jì)精度的先驗(yàn)信息,并把這些先驗(yàn)信息合理地融入當(dāng)前已有的反演成像算法中。對實(shí)際應(yīng)用而言,定量地震波成像對于成像處理各環(huán)節(jié)是否做到了保真的要求很高,因此質(zhì)量控制是十分重要的。一項(xiàng)好的方法技術(shù)要落地應(yīng)用,不僅僅是靠理念正確和方法技術(shù)的有效,各環(huán)節(jié)的質(zhì)控往往起到?jīng)Q定成敗的作用。

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

猜你喜歡
子波波阻抗波數(shù)
聲場波數(shù)積分截?cái)嗖〝?shù)自適應(yīng)選取方法
一種基于SOM神經(jīng)網(wǎng)絡(luò)中藥材分類識(shí)別系統(tǒng)
電子測試(2022年16期)2022-10-17 09:32:26
一類非線性動(dòng)力系統(tǒng)的孤立子波解
波阻抗技術(shù)在煤礦三維地震勘探中的應(yīng)用
波阻抗使用單位規(guī)范問題探究
地震反演子波選擇策略研究
波阻抗反演技術(shù)與砂體理論模型的對比
重磁異常解釋的歸一化局部波數(shù)法
基于聲場波數(shù)譜特征的深度估計(jì)方法
基于倒雙譜的地震子波估計(jì)方法
扎兰屯市| 嵊泗县| 庆安县| 上犹县| 新密市| 宝清县| 叙永县| 磴口县| 阳江市| 拉萨市| 绥棱县| 濮阳县| 新乡县| 靖西县| 大足县| 门头沟区| 怀来县| 昆山市| SHOW| 荥经县| 仁布县| 芦山县| 池州市| 涿州市| 那曲县| 沂水县| 平顶山市| 长乐市| 四子王旗| 竹北市| 陆河县| 苏尼特右旗| 大同市| 固始县| 藁城市| 汶上县| 八宿县| 留坝县| 鹤岗市| 新民市| 吉木萨尔县|