劉玉柱,楊積忠
(1.同濟大學(xué)海洋地質(zhì)國家重點實驗室,上海200092;2.同濟大學(xué)海洋與地球科學(xué)學(xué)院,上海200092)
在山地或復(fù)雜表層探區(qū)地震勘探中,精細的表層速度結(jié)構(gòu)對于靜校正、基準(zhǔn)面延拓、起伏地表偏移成像,甚至全波形反演等地震數(shù)據(jù)處理至關(guān)重要。在傳統(tǒng)的勘探地震學(xué)方法中,主要采用基于射線理論的地震層析成像技術(shù),利用初至波或折射波旅行時反演近地表宏觀速度結(jié)構(gòu)[1-2],然而這樣獲得的速度結(jié)構(gòu)已經(jīng)無法滿足目前后續(xù)地震數(shù)據(jù)處理流程對高精度近地表速度結(jié)構(gòu)的要求。為此,很多學(xué)者從層析成像理論上探索提高表層速度建模精度的方法,如劉玉柱等[3]放棄傳統(tǒng)的射線理論,基于有限頻理論提出了初至波菲涅爾體地震走時層析成像方法。利用該方法可以有效提高表層建模精度,但計算效率顯著降低。Sheng等[4]基于波動理論提出了初至波波形反演方法。該方法在得到精細的表層速度結(jié)構(gòu)的同時,可以降低波形反演的多解性與不確定性,但由于需要理論合成地震波場,其計算效率更低,且反演結(jié)果受地震資料信噪比、地震子波未知等眾多實際因素的影響。
利用地面觀測數(shù)據(jù)反演高精度的近地表速度結(jié)構(gòu)異常困難,因為除以上層析成像理論本身的局限外,還有很多因素影響層析效果[5],反演策略就是其中一個重要的影響因素。為此,很多學(xué)者對層析反演策略進行了研究。如為了同時反演大尺度與小尺度速度結(jié)構(gòu),Zhou[6]提出了基于多重網(wǎng)格的多尺度地震層析成像方法;為了有效利用模型先驗信息,劉玉柱等[7]提出了針對不同先驗信息的正則化處理方法;為了解決復(fù)雜起伏地表的層析成像問題,Zhou等[8]在多重網(wǎng)格層析成像策略的基礎(chǔ)上進一步提出了約束可形變層析成像方法;為了利用地震數(shù)據(jù)的多震相信息,Zhou[9]利用地表VSP初至數(shù)據(jù)反演表層速度結(jié)構(gòu),Li等[10]提出了聯(lián)合初至與淺層反射的聯(lián)合層析成像方法。然而,這些研究針對的是如何利用模型先驗信息或者如何利用更多的數(shù)據(jù)反演表層速度結(jié)構(gòu),關(guān)于如何有效利用數(shù)據(jù)中包含的先驗信息的研究很少。
根據(jù)不同的反演目標(biāo)制定不同的反演策略,在提高反演精度的同時可以達到事半功倍的效果。我們根據(jù)不同深度處反演精度與不同偏移距初至數(shù)據(jù)的實驗關(guān)系,提出有效利用多偏移距數(shù)據(jù)的偏移距加權(quán)地震層析成像方法,以提高表層速度結(jié)構(gòu)的反演精度。
圖1為二維復(fù)雜起伏地表理論模型。該模型包含4001×151個網(wǎng)格,網(wǎng)格大小為10m×10m,速度變化范圍為800~5500m/s?;谠撃P瓦M行正演模擬,第1個激發(fā)點位于水平方向5000m處,激發(fā)點水平間隔40m,中間激發(fā)兩邊接收,接收點水平間距為20m,最大偏移距4000m,最小偏移距為0。根據(jù)該觀測系統(tǒng)進行最小走時射線追蹤,利用彈性波方程合成了640炮記錄。
圖1 二維復(fù)雜起伏地表理論模型(白色為空氣層)
圖2為模擬數(shù)據(jù)不同偏移距范圍內(nèi)的射線覆蓋密度圖。從圖2可以看出,小偏移距數(shù)據(jù)對應(yīng)的射線路徑主要集中在淺層,而且分布比較均勻。隨著偏移距范圍的增大,雖然射線路徑對淺層與深部都有覆蓋,但對深部的覆蓋次數(shù)明顯高于淺層,這說明深部對大偏移距數(shù)據(jù)的影響更大。
為了進行更加細致的對比,圖3給出了理論模型中x=20km處所有深度上不同偏移距范圍數(shù)據(jù)對應(yīng)的射線覆蓋密度圖。圖3更加清晰地表明:
小偏移距數(shù)據(jù)對應(yīng)的射線路徑主要集中在淺層,且在淺層的覆蓋次數(shù)總體上要大于大偏移距數(shù)據(jù)在淺層的覆蓋次數(shù);大偏移距數(shù)據(jù)對應(yīng)的射線路徑主要集中在深層,且總體上在深層具有絕對的覆蓋,即小偏移距數(shù)據(jù)的射線路徑基本上不經(jīng)過深層。
根據(jù)射線層析成像成功應(yīng)用的假設(shè)前提,即只有射線路徑經(jīng)過的空間區(qū)域?qū)?yīng)的模型參數(shù)才能夠被反演出來,而且反演的精度與射線的覆蓋密度成正比,我們可以得出如下結(jié)論:①對于基于地表觀測系統(tǒng)的地震數(shù)據(jù)初至層析成像而言,淺層的精細速度結(jié)構(gòu)隱藏在所有偏移距數(shù)據(jù)當(dāng)中,且主要隱藏在小偏移距的初至數(shù)據(jù)當(dāng)中,深層的速度結(jié)構(gòu)則完全隱藏在大偏移距的初至數(shù)據(jù)當(dāng)中;②利用大偏移距的初至數(shù)據(jù)可以同時反演得到淺層及深層的模型參數(shù),但淺層的反演精度比較低;③欲反演淺層的精細速度結(jié)構(gòu),在利用大偏移距初至數(shù)據(jù)的同時,需重點利用小偏移距初至數(shù)據(jù),欲反演深層速度結(jié)構(gòu)則需要重點利用大偏移距初至數(shù)據(jù)。
圖2 基于理論模型的不同偏移距范圍數(shù)據(jù)對應(yīng)的射線覆蓋密度a -500~500m; b -1000~1000m; c -2000~2000m; d -4000~4000m
圖3 理論模型中x=20km處所有深度上不同偏移距范圍數(shù)據(jù)對應(yīng)的射線覆蓋密度堆疊顯示
在線性假設(shè)前提下,基于射線理論的初至波走時層析成像可簡化為求解大規(guī)模線性方程組
(1)
其中,L為m×n維矩陣,矩陣元素為lkj,代表第k條射線在第j個模型參數(shù)單元內(nèi)的長度;Δs為n維列向量,代表模型慢度參數(shù)修正量;Δt為m維列向量,代表觀測初至波走時相對于理論計算初至波走時的殘差[11]。層析方程中傳播矩陣L的構(gòu)建是反演問題的關(guān)鍵。(1)式所對應(yīng)的阻尼最小二乘解如下[7]:
(2)
為了重點利用小偏移距數(shù)據(jù)反演表層精細速度結(jié)構(gòu),基于第1章數(shù)值實驗的結(jié)果,我們提出多偏移距(multi-offset range,MOR)層析反演策略。首先根據(jù)大偏移距的初至信息進行層析反演,獲得所有深度的宏觀速度分布;然后將其作為初始模型,再以較小偏移距的初至信息進行層析反演,以提高淺層的反演精度;如此反復(fù),直到反演精度滿足要求為止。為便于區(qū)分,我們將傳統(tǒng)的基于(1)式和(2)式的單次層析成像反演策略稱為單偏移距(single-offset range,SOR)層析策略。
圖4是理論模型(圖1)數(shù)據(jù)的MOR層析反演結(jié)果,圖5是某實際二維測線資料的MOR層析反演結(jié)果。不難看出,無論是模型數(shù)據(jù)還是實際資料的反演處理,每縮小一次偏移距范圍,表層的速度結(jié)構(gòu)就更加精細一些,最終使表層的低速結(jié)構(gòu)被成功反演了出來,這是傳統(tǒng)的SOR初至波層析成像方法所無法實現(xiàn)的。理論模型與實際資料在證實了MOR層析反演策略有效性的同時,也驗證了第1章實驗所得結(jié)論的正確性。
圖4 理論模型數(shù)據(jù)MOR層析反演結(jié)果a 利用-2000~2000m偏移距范圍內(nèi)的初至走時; b 以圖4a為初始模型,利用-1000~1000m偏移距范圍內(nèi)的初至走時; c 以圖4b為初始模型,利用-500~500m偏移距范圍內(nèi)的初至走時
圖5 某實際二維測線資料MOR層析反演結(jié)果a 利用-2000~2000m偏移距范圍內(nèi)的初至走時; b 以圖5a為初始模型,利用-1000~1000m偏移距范圍內(nèi)的初至走時; c 以圖5b為初始模型,利用-500~500m偏移距范圍內(nèi)的初至走時; d 以圖5c為初始模型,利用-250~250m偏移距范圍內(nèi)的初至走時
盡管MOR層析反演策略可以有效利用小偏移距地震數(shù)據(jù)反演得到表層的精細速度結(jié)構(gòu),但該方法多次嵌套反演的過程繁瑣,不易操作。
根據(jù)Tarantola和Valette[12]對非線性反演問題目標(biāo)函數(shù)的定義可以得到初至波走時層析方程(1)滿足方程(3)的最優(yōu)解表達式
(3)
(4)
(5)
(6)
式中:hmax為最大偏移距;hsg表示當(dāng)前偏移距;d為最大權(quán)系數(shù)(即零偏移距數(shù)據(jù)對應(yīng)的權(quán)系數(shù)),是一個比較大的正數(shù)。公式(6)表明,最大偏移距的權(quán)系數(shù)為d/(d+1),d較大時其接近于1,而零偏移距的權(quán)系數(shù)為d。
相比于MOR層析反演,偏移距加權(quán)地震層析成像反演中的權(quán)重系數(shù)是偏移距的光滑函數(shù)而非MOR中的跳躍函數(shù),因此該方法可以實現(xiàn)對所有偏移距數(shù)據(jù)的區(qū)別利用。同時,該方法只需一次SOR反演即可實現(xiàn)MOR層析中的多次嵌套反演,大大簡化了反演流程。
為了驗證偏移距加權(quán)地震層析成像方法的有效性,本文基于上述理論模型與實際資料再次進行了層析成像反演實驗。在實驗處理中,取最大權(quán)系數(shù)d=10。圖6和圖7給出了分別基于上述理論模型數(shù)據(jù)和實際資料進行偏移距加權(quán)地震層析成像反演的結(jié)果(其中圖6為表層的放大顯示)。
由圖6可以看出,理論模型表層的低速結(jié)構(gòu)基本上被反演了出來;相比于2000m偏移距范圍內(nèi)的SOR反演結(jié)果(圖4a),圖6所示結(jié)果具有明顯的優(yōu)勢,反演精度與分辨率更高;從總體上看,圖6所示結(jié)果稍遜于MOR反演結(jié)果(圖4c),但有些區(qū)域反演結(jié)果(圖中方框標(biāo)識的區(qū)域)則更加接近理論模型(圖1)。
圖6 理論模型數(shù)據(jù)偏移距加權(quán)地震層析成像反演結(jié)果
圖7所示實際資料反演結(jié)果與MOR層析反演結(jié)果(圖5d)相當(dāng),表層的低速精細結(jié)構(gòu)基本都反演了出來。在圖7中矩形與橢圓標(biāo)識處,偏移距加權(quán)層析反演結(jié)果比MOR層析反演結(jié)果更加精細。兩者一個比較大的區(qū)別在于,圖7中箭頭所示處是一個比較大的低速體,而MOR反演結(jié)果顯示該處(圖5d箭頭所示)是一個較小的低速體。這是由于這兩種方法都是針對表層的速度反演,在該處的反演結(jié)果具有較大的不確定性。
圖7 某實際二維測線資料偏移距加權(quán)地震層析成像反演結(jié)果
經(jīng)過多個理論模型與實際資料的測試發(fā)現(xiàn),最大權(quán)系數(shù)d的選取需要注意:公式(3)是利用數(shù)據(jù)協(xié)方差矩陣中的誤差實現(xiàn)對數(shù)據(jù)的加權(quán)利用,而在本文中為了反演表層速度結(jié)構(gòu),將公式(3)的數(shù)學(xué)含義轉(zhuǎn)化為物理含義,通過偏移距加權(quán)實現(xiàn)對小偏移距數(shù)據(jù)的重點利用。因此,對于實際資料而言,在無法確定觀測數(shù)據(jù)誤差、但可以確定誤差與偏移距是正相關(guān)的條件下,d應(yīng)該選取較大的數(shù);在不考慮觀測誤差的情況下,d可以選取一個相對較小的數(shù)。
基于初至層析反演中不同偏移距數(shù)據(jù)與不同深度反演精度的實驗認識,提出了MOR層析反演的策略,并通過理論模型與實際資料對其有效性進行了驗證。為了簡化MOR層析反演的迭代步驟,根據(jù)數(shù)據(jù)協(xié)方差矩陣在廣義反演目標(biāo)函數(shù)中的物理意義,進一步給出了新的數(shù)據(jù)協(xié)方差矩陣表達形式,進而提出了偏移距加權(quán)地震層析成像方法。通過理論模型數(shù)據(jù)與實際資料處理實驗證實了該方法的有效性。
在偏移距加權(quán)地震層析成像方法中,數(shù)據(jù)的權(quán)重系數(shù)是偏移距的光滑函數(shù)而非MOR中的跳躍函數(shù),因此可以實現(xiàn)對所有偏移距數(shù)據(jù)的區(qū)別利用。同時,該方法只需一次反演即可實現(xiàn)MOR中的多次嵌套反演,大大簡化了反演流程。但由于偏移距加權(quán)函數(shù)不是從理論上導(dǎo)出,而是基于經(jīng)驗認識給出的等效公式,因此加權(quán)函數(shù)本身及最大權(quán)系數(shù)d的選取值得進一步研究。
參 考 文 獻
[1] Zhu X H,Sixta D P,Angstman B G.Tomostatics:turning-ray tomography+static corrections [J].The Leading Edge,1992,11(12):15-23
[2] 李錄明,羅省賢,趙波.初至波表層模型層析反演[J].石油地球物理勘探,2000,35(5):559-564
Li L M,Luo X X,Zhao B.Tomographic inversion of first break in surface model [J].Oil Geophysical Prospecting,2000,35(5):559-564
[3] 劉玉柱,董良國,王毓瑋,等.初至波菲涅爾體地震層析成像[J].地球物理學(xué)報,2009,52(9):2310-2320
Liu Y Z,Dong L G,Wang Y W,et al.Fresnel volume tomography based on the first arrival of the seismic wave [J].Chinese Journal of Geophysics,2009,52(9):2310-2320
[4] Sheng J M,Alan L,Maike B,et al.Early arrival waveform tomography on near-surface refraction data[J].Geophysics,2006,71(4):47-57
[5] 劉玉柱,董良國.近地表速度結(jié)構(gòu)初至波層析影響因素分析[J].石油地球物理勘探,2007,42(5):544-553
Liu Y Z,Dong L G.Analysis of influence factors of first-breaks tomography [J].Oil Geophysical Prospecting,2007,42(5):544-553
[6] Zhou H W.Multiscale traveltime tomography [J].Geophysics,2003,68(5):1639-1649
[7] 劉玉柱,董良國,夏建軍.初至波走時層析中的正則化方法[J].石油地球物理勘探,2007,42(6):682-685,698
Liu Y Z,Dong L G,Xia J J.Normalized approach in tomographic imaging of first breaks travel time [J].Oil Geophysical Prospecting,2007,42(6):682-685,698
[8] Zhou H W,Li P M,Yan Z H,et al.Constrained deformable layer tomostatics [J].Geophysics,2009,74(6):WCB35-WCB46
[9] Zhou H W.First-break vertical seismic profiling tomography for Vinton Salt Dome [J].Geophysics,2006,71(3):U29-U36
[10] Li P M,Yan Z H,Guo M J,et al.2D deformable-layer tomostatics with the joint use of first breaks and shallow reflections [J].Expanded Abstracts of 79th
Annual Internat SEG Mtg,2009,1335-1339
[11] Stork C,Clayton R W.Linear aspects of tomographic velocity analysis [J].Geophysics,1991,56(4):483-495
[12] Tarantola A,Valette B.Generalized non-linear inverse problems solved using the least-squares criterion [J].Reviews of Geophysics and Space Physics,1982,20(2):219-232