地層格架混合總則化及其應(yīng)用

0基于層析的地質(zhì)結(jié)構(gòu)模型求解應(yīng)用于偏移成像的宏觀速度模型構(gòu)建是反射地震圖像處理的重要步驟。基于傳統(tǒng)的成像環(huán)境分析反演速度建模，已進入行業(yè)標(biāo)準(zhǔn)偏移速度建模流程。與傳統(tǒng)的成像環(huán)境分析方法不同，三維anova方法是數(shù)據(jù)域分析的一種方法。與以往基于時間的反射分析算法相比，該算法充分利用地震波場中所有運動信息。將本地相位軸在炮道收集和檢測點收集的波束參數(shù)（以下簡稱p參數(shù)）的水平分量、炮點坐標(biāo)和檢測點坐標(biāo)納入分析反演的數(shù)據(jù)空間，并重新排列模型空間和數(shù)據(jù)空間。三維anova反演是唯一一種抗彎折射方法，可以同時反演速度、反射點位置和局部地層傾角。對于工業(yè)界用戶來說,最經(jīng)典的獲得地質(zhì)一致性塊狀模型的方式當(dāng)屬層剝離方法.該方法從20世紀(jì)90年代初即進入商業(yè)軟件,至今仍為工業(yè)界廣泛應(yīng)用.但是層剝離方法的問題是淺層的速度和界面位置一旦確定就無法再改動,誤差會從淺層慢慢積累至深層,導(dǎo)致深層的速度建模產(chǎn)生大的偏差.因此有學(xué)者開始探索是否可以在射線類層析方法中同時實現(xiàn)界面和速度的更新.近年來在射線類層析的研究進展方面,有兩類獲得地質(zhì)一致性塊狀模型的方法值得關(guān)注.第一種方法是使用混合正則化項,比如在地質(zhì)界面處使用L1范數(shù)或其他稀疏類范數(shù)使得反演結(jié)果保留模型的突變特征,而在其他地方則使用L2范數(shù)或其他光滑類范數(shù)保持反演結(jié)果的平滑(基于射線類層析方法獲得地質(zhì)一致性塊狀結(jié)構(gòu)速度的目標(biāo)非常明確:反演得到的速度模型應(yīng)在保證地質(zhì)塊內(nèi)足夠光滑的同時盡可能維持邊界的突變.對于立體層析而言,這個目標(biāo)可以通過更簡潔的方法來實現(xiàn).本文提出一種實用的模型正則化方法,該方法的核心思想是將封閉后的地層格架信息融入到立體層析Fréchet導(dǎo)數(shù)矩陣中,在通過融入地層格架信息簡化了立體層析Fréchet導(dǎo)數(shù)矩陣之后,使得更新后的速度模型自動呈現(xiàn)出與地質(zhì)構(gòu)造較為一致的的層狀或者塊狀結(jié)果.其實現(xiàn)過程如下:(1)通過非均勻B樣條方法對更新后的反射點位置進行自動擬合得到更新后地質(zhì)層位信息;(2)利用地質(zhì)塊體之間的拓撲關(guān)系對地質(zhì)層位信息實施自動封閉得到封閉后的地層格架信息;(3)利用封閉后的地層格架信息對Fréchet導(dǎo)數(shù)矩陣進行分塊,每塊內(nèi)的速度更新量以及速度梯度更新量均為常數(shù),由于一個模型內(nèi)的地質(zhì)塊數(shù)目一般不會很多,這樣處理的結(jié)果使得模型空間中的速度分量數(shù)目被大幅度壓縮,并且在更新過程中模型自動呈現(xiàn)出層狀或者塊狀地質(zhì)特征.設(shè)計上述實現(xiàn)策略的原因是:反射點位置本身即屬于立體層析模型空間的一部分,反射點位置的更新在立體層析反演的迭代過程中是可以自動得到的.相比在每一輪反演后實施偏移成像再進行人工拾取,這是一種不僅效率更高、而且更為合理的方法.同時相比于本文第1節(jié)首先介紹傳統(tǒng)立體層析方法的數(shù)據(jù)空間和模型空間以及本文提出的地層格架約束下立體層析的數(shù)據(jù)空間和模型空間;第2節(jié)給出了地層格架正則化下的立體層析Fréchet導(dǎo)數(shù)矩陣的具體算法;第3節(jié)給出了基于散亂反射點的非均勻B樣條擬合算法與地層格架自動封閉的實現(xiàn)過程;第4節(jié)給出一個實用的工作流程;第5節(jié)給出了理論與實際數(shù)據(jù)應(yīng)用算例.作者期望這種特殊的正則化方式獲得的立體層析反演結(jié)果將可以用于逆時偏移處理或作為全波形反演的高質(zhì)量初始模型.1有限層框架下的二維三維層析反演1.1維立體層析freget偏導(dǎo)數(shù)矩陣公式首先介紹常規(guī)二維立體層析的模型空間、數(shù)據(jù)空間及Fréchet導(dǎo)數(shù)矩陣矩陣的建立過程.圖1a顯示了一根從炮點S出發(fā)、到地下反射點X反射、回到地表R的射線.圖1b顯示了射線中心坐標(biāo)系中射線起點和端點擾動的幾何關(guān)系.從透射的角度,不妨將其理解為從X出發(fā)、分別以入射角θ任何一種層析反演方法都必須建立數(shù)據(jù)空間擾動與模型空間擾動的線性關(guān)系,即:這樣利用觀測到的數(shù)據(jù)殘差就可以計算出模型空間的擾動量,達到更新初始模型的目的.其中關(guān)鍵是建立方程(1)中所示的F矩陣,即Fréchet偏導(dǎo)數(shù)矩陣公式.Fréchet偏導(dǎo)數(shù)矩陣的物理意義是數(shù)據(jù)空間任一分量關(guān)于模型空間任一分量的敏感度.在常規(guī)走時層析中,走時t關(guān)于速度V的偏導(dǎo)數(shù)其實就是射線經(jīng)過每一個速度單元的弧長.(2)式展示了二維常規(guī)立體層析Fréchet偏導(dǎo)數(shù)矩陣,除了第一行走時關(guān)于模型空間的偏導(dǎo)數(shù)根據(jù)走時積分方程即可得到外.其余元素是都通過射線擾動理論(Farraand1.2層位優(yōu)化設(shè)計在傳統(tǒng)立體層析中速度模型一般由B樣條基函數(shù)的權(quán)系數(shù)來描述,這種表達可以達到壓縮立體層析Fréchet導(dǎo)數(shù)矩陣規(guī)模的目的.然而,B樣條描述方式的一個副作用是會使得反演結(jié)果過于平滑.與此同時,傳統(tǒng)的正則化技術(shù)一般都建議加強解的平滑性通過一個概念實驗即可清楚闡述地層格架約束立體層析的核心思想.圖2顯示了一個6層背斜模型.如果按照均勻B樣條基函數(shù)表達這個模型,在節(jié)點縱橫向間距都為400m(一般取200m或400m)時也需要480個B樣條系數(shù),也就是說模型空間中關(guān)于速度的未知數(shù)有480個.如果說用常規(guī)的規(guī)則網(wǎng)格表達圖2所示的速度模型,在網(wǎng)格縱橫向間隔都設(shè)為100m時,會有120×48=5760個網(wǎng)格,也就是說模型空間中關(guān)于速度的未知數(shù)將有5760個.然而,如果我們假定每一塊內(nèi)的速度都應(yīng)該是常數(shù)時,其模型空間內(nèi)的速度未知數(shù)其實應(yīng)該是6個.基于一個由6個地質(zhì)界面封閉得到的地層格架實施上述規(guī)則化就相當(dāng)于將速度未知量的數(shù)目從一個很大的規(guī)模壓縮為6個.同時,由于每一塊內(nèi)的速度更新量各有差異,因此更新后的模型在地質(zhì)層位兩側(cè)必然會呈現(xiàn)出強反差.這正是地層格架約束立體層析所想要得到的地質(zhì)一致性塊狀特征.一個關(guān)鍵問題是:地層格架如何在反演中的到正確的更新?注意地層格架是由許多地質(zhì)層位封閉成的,而這些地質(zhì)層位是由許多反射點構(gòu)建的.而這些反射點的位置本身就是立體層析模型空間的一部分,在每一次迭代中它們都和速度、散射角等模型信息一起被更新.因此完全可以利用更新后的反射點位置來擬合出更新后的地質(zhì)層位.這里我們使用非均勻B樣條(NUBS)(2維中心坐標(biāo)系下地表格架內(nèi)射線擾動的一階關(guān)系那么二維地層格架約束立體層析的Fréchet導(dǎo)數(shù)如何計算?注意二維地層格架約束立體層析的數(shù)據(jù)空間與傳統(tǒng)二維立體層析并無差別,依然可以描述為d=(S炮點一側(cè)的數(shù)據(jù)空間對于檢波點一側(cè)的散射角顯然沒有敏感度,反之亦然.同時根據(jù)費馬原理,走時T在速度、炮點與接收點坐標(biāo)都固定時已經(jīng)是最小走時,因此散射角對走時亦不存在一階擾動,有:根據(jù)Farra和Madariaga(1987),二維中心坐標(biāo)系下立體層析所用的線性擾動關(guān)系如(4)、(5)式所示.(4)式等號左邊為地表觀測的坐標(biāo)擾動(Δq(4)式中的Δq將(7)式代入(4)式容易推出:方括號中后兩項為顯然,在射線到達地表后,需要將擾動量從中心坐標(biāo)系換算到地表水平觀測面.因此需要建立射線中心坐標(biāo)系下的位移擾動Δq將(6)、(7)、(8)、(9)、(10)式代入到(4)式中,可得出二維中心坐標(biāo)系下地表觀測數(shù)據(jù)信息與初始射線擾動擾動信息的一階關(guān)系為如前所述,當(dāng)速度模型被地層格架信息分塊之后,速度更新量在一個地質(zhì)塊內(nèi)將是一個常數(shù),因此(12)式中括號里的第一項,即速度更新量關(guān)于空間坐標(biāo)的導(dǎo)數(shù)就應(yīng)該是零;其次一個地質(zhì)塊內(nèi)的速度定義為根據(jù)鏈?zhǔn)椒▌t,二維中心坐標(biāo)系下坐標(biāo)擾動Δq、慢度擾動Δp關(guān)于v最終我們得到:類似地,將鏈?zhǔn)椒▌t應(yīng)用于(5)式,就可以得到走時關(guān)于v至此,在引入地層格架信息約束后,走時t,地表觀測坐標(biāo)x,地表觀測水平慢度矢量p3基于非均勻b樣條的地質(zhì)建模對地層格架約束立體層析而言,在反演實施之前對散亂反射點的擬合和自動封閉至關(guān)重要.這里我們采用了非均勻B樣條擬合算法將散亂反射點擬合為連續(xù)的地質(zhì)層位.相比均勻B樣條擬合,該算法能夠高效擬合多值函數(shù),其穩(wěn)健性在其他行業(yè)中已經(jīng)得到了充分體現(xiàn),證明完全能夠適應(yīng)二維、三維復(fù)雜地質(zhì)建模的需求.這里對其實現(xiàn)原理做一簡介,更多細節(jié)請讀者參閱(Piegland3.1未知控制點求解B樣條基函數(shù)的定義如下:令U={u而對于p次B樣條曲線的定義則為對于數(shù)據(jù)擬合,需要預(yù)先計算好數(shù)據(jù)點的參數(shù)值和節(jié)點矢量,然后,建立并求解線性最小二乘問題來求解未知控制點.假定p≥1,n≥p,并且給定數(shù)據(jù)點Q其滿足條件:Q求取擬合曲線主要分為以下幾個步驟:(1)預(yù)先通過弦長參數(shù)化得到數(shù)據(jù)點Q由于節(jié)點的分布應(yīng)該反映節(jié)點矢量的分布,假定c是一個正實數(shù),用i=int(c)表示小于或等與d的最大整數(shù),則總共需要n+p+2個節(jié)點,因此有n-p個內(nèi)節(jié)點和n-p-1個內(nèi)部節(jié)點區(qū)間,令c=(m+1)/(n-p+1),然后按公式(20)定義內(nèi)節(jié)點,即:(2)建立目標(biāo)函數(shù),公式為3.2層位組合格架封閉第3.1節(jié)中的非均勻B樣條擬合僅僅是針對某一個地質(zhì)層位的離散反射點完成了連續(xù)地質(zhì)層位的擬合.但是要實施地層格架的自動封閉必須考慮地質(zhì)層位之間的相交或者令其自然延伸到剖面的邊界處,這些在實施自動封閉時都需要非常仔細的考慮.此外,對于含有多個地質(zhì)塊體的地震剖面,地層的封閉不但需要考慮塊體頂?shù)酌嬷g的關(guān)系,還要考慮塊體與地層以及地層與地層之間的關(guān)系.這要求離散的反射點數(shù)據(jù)必須含有地質(zhì)層位信息,即在反演之前,在自動拾取主要地質(zhì)層位的立體層析數(shù)據(jù)空間時就需要加入層位編號.實施地層格架封閉分為以下5個步驟:(1)獲取層位編號信息,獲取塊體與地層、地層與地層之間的幾何拓撲關(guān)系,即確定塊體與層位之間、層位與層位之間的拓撲關(guān)系(比如一個交點由哪些層位組成,一個塊體由哪些層位包圍).(2)通過反射點信息擬合各個層位:對每個層位的數(shù)據(jù)點根據(jù)橫坐標(biāo)進行排序后,然后對同一個層位的數(shù)據(jù)進行NUBS擬合,利用結(jié)構(gòu)張量算法提取擬合后界面的切向量和法向量信息(Wuand(3)處理界面與邊界之間的關(guān)系:在確定層位之間的關(guān)系時,統(tǒng)計每個層位的交點個數(shù)(只有0、1、2三種情況),當(dāng)只有一個交點時,判斷哪個端點需要延拓到邊界,當(dāng)沒有交點時,地層兩端根據(jù)當(dāng)前的界面端點處的斜率延拓到邊界,斜率信息是通過結(jié)構(gòu)張量算法獲得的.(4)處理界面之間的拓撲關(guān)系:對于每個交點,確定哪些界面相交于該點,然后通過一些處理使得這些層位相交于一點,處理過程中同樣需要結(jié)構(gòu)張量獲取的界面的切向量和法向量信息.(5)根據(jù)塊體與層位之間的拓撲關(guān)系求取塊體頂?shù)酌娴纳疃茸鴺?biāo),每一個塊體內(nèi)部一個橫坐標(biāo)的頂?shù)酌嫔疃仁菍崿F(xiàn)分塊求取Fréchet導(dǎo)數(shù)的重要參數(shù).圖4顯示了基于離散的反射點實施地層格架封閉的全過程.圖4a顯示了塊體與界面之間、界面與界面之間的拓撲關(guān)系;圖4b顯示了如何通過反射點信息擬合各個反射界面;圖4c顯示了將界面延伸到邊界后的情況;圖4d顯示了最為重要也是相對復(fù)雜的一步:處理好界面與界面之間的交點關(guān)系.4基于離散反射點的連續(xù)界面擬合到目前為止我們已經(jīng)實現(xiàn)了:(1)推導(dǎo)了二維地層格架約束立體層析的Fréchet導(dǎo)數(shù);(2)利用NUBS插值實現(xiàn)基于離散反射點的連續(xù)界面擬合;(3)利用初始成像界面提供的地質(zhì)塊與界面的拓撲結(jié)構(gòu)實現(xiàn)地層格架的自動封閉.根據(jù)上述這些算法可以得到一個實用的、適用于二維實際數(shù)據(jù)的工作流程.該流程實現(xiàn)過程如下:(1)對輸入的疊前地震數(shù)據(jù)利用結(jié)構(gòu)張量算法提取高密度立體層析數(shù)據(jù)空間[S5理論數(shù)據(jù)和真實數(shù)據(jù)的計算5.1初始速度模型的建立和反演的驗證圖6顯示這個工作流程在圖2所示的六層背斜模型理論數(shù)據(jù)上是如何工作.基于圖2顯示的模型,共正演了6000個射線對并得到對應(yīng)的立體層析數(shù)據(jù)空間.首先運行一個傳統(tǒng)立體層析得到初始反射點位置和初始速度模型(如圖6a所示).圖6b顯示基于圖6a中獲得的初始反射點位置實施NUBS擬合可以獲得6個連續(xù)的反射界面.圖6c顯示根據(jù)圖6b所示反射界面實施地層格架封閉后的結(jié)果.地層格架約束立體層析的初始模型是將圖6c所示的地層格架內(nèi)充填一個初始速度后得到的.初始充填速度的獲得方式和圖5的建議流程完全一致,即采取在圖6a所示的初始模型內(nèi)每一個地質(zhì)塊內(nèi)做一個簡單的算術(shù)平均后得到的.圖6d顯示了地層格架約束二維立體層析在這個模型上的最終反演結(jié)果.我們可以看出圖6d所示的反演結(jié)果和圖2所示的真實模型已經(jīng)非常接近,證明這個流程在圖2所示的層狀模型上是有效的.5.2基于深度成像的地層格架和控制律設(shè)計的優(yōu)化模型圖7顯示這個工作流程在另一個更為復(fù)雜的理論數(shù)據(jù)上如何工作.基于圖7a顯示的模型,正演了一條2D理論測線,共激發(fā)201炮,炮間距為40m,每炮401道接收,道間距為10m,最大偏移距2km.在自動搜索到5700個數(shù)據(jù)點后,首先運行一個傳統(tǒng)立體層析得到初始反射點位置和初始速度模型(如圖7b).基于圖7b顯示的反射點位置通過NUBS擬合算法可以實現(xiàn)地層格架自動封閉(如圖7d所示).圖7e顯示了基于圖7b顯示的初始模型和圖7d所示的封閉后地層格架在每一塊內(nèi)做一個簡單的平均后得到了用于二維地層格架約束立體層析的初始模型.在圖7e基礎(chǔ)上,我們運行了地層格架約束立體層析.圖7f、g分別顯示地層格架約束立體層析第1次和第10次迭代更新的速度模型.我們可以看出圖7g所示的最終反演結(jié)果和圖7a所示的真實模型已經(jīng)非常接近.圖7h為目標(biāo)函數(shù)的下降曲線,可以看出最終的二范數(shù)殘差已經(jīng)非常小.圖8a、b顯示了用圖7b所示傳統(tǒng)立體層析模型實施偏移成像得到的共成像點道集(CIG)以及用圖7g所示地層格架約束立體層析反演模型實施偏移成像獲得的CIG之間的對比.圖8c、d顯示了用圖7b所示傳統(tǒng)立體層析模型實施偏移成像得到的成像剖面以及用圖7g所示地層格架約束立體層析反演模型實施偏移成像獲得的成像剖面的對比.可以看出,本文方法獲得的深度成像剖面的聚焦程度和成像品質(zhì)相比傳統(tǒng)立體層析有明顯提高.5.3基于數(shù)據(jù)信噪比的地層格架約束立體層析二維實際數(shù)據(jù)是2011年在南海深水區(qū)采集獲得,共選取了該測線中950炮,炮間距為50m,道間距為25m,最大偏移距為8275m.考慮到實際模型的復(fù)雜性,每一層內(nèi)的速度設(shè)置為隨縱、橫向線性變化.圖9a顯示了傳統(tǒng)立體層析得到的最終反演模型,其中反演得到的傾角條用藍色線畫在速度模型之上.應(yīng)用同樣的工作流程,在NUBS擬合算法和結(jié)構(gòu)張量算法的幫助下,我們得到了一個封閉的地質(zhì)格架.由于該數(shù)據(jù)信噪比的原因,僅僅基于了幾個關(guān)鍵層位建立了一個比較粗放的地層格架,沒有針對更精細的反射點信息實施地層格架封閉.盡管如此,地層格架約束立體層析對于深部的速度結(jié)構(gòu)反演依然有明顯改善.圖9b、c顯示地層格架約束立體層析的第1次和第7次迭代結(jié)果.圖8d顯示迭代過程中目標(biāo)函數(shù)的下降曲線.由于實際數(shù)據(jù)存在信噪比以及數(shù)據(jù)點提取的密度和精度問題,因此其目標(biāo)函數(shù)的下降幅度不可能像理論數(shù)據(jù)的目標(biāo)函數(shù)那樣,可以下降到接近于零的水平.但是在引入地層格架之后,對于測線中部5000～8000m深度的內(nèi)幕的成像依然有明顯改進.圖9e、f顯示了用圖9a所示傳統(tǒng)立體層析模型實施偏移成像得到的共成像點道集(CIG)以及用圖9c所示地層格架約束立體層析反演模型實施偏移成像獲得的CIG之間的對比.圖9g