基于正態分布權函數的VSP-CDP疊加成像

楊飛龍 李輝峰 孫 輝 張 雪 羅 浩 趙 馳

(①西安石油大學地球科學與工程學院，陜西西安710065； ②陜西省油氣成藏地質學重點實驗室，陜西西安710065；③山東科技大學地球科學與工程學院山東省沉積成礦作用與沉積礦產重點實驗室，山東青島 266590；④西南交通大學地球科學與環境工程學院，四川成都 611756)

0 引言

隨著油氣勘探的深入，面對規模小、結構復雜、物性變化劇烈的地層圈閉和巖性圈閉，對地震勘探分辨率的要求越來越高。垂直地震剖面法勘探(VSP)是一種井中地震觀測技術，主要研究地球介質內部波場特征、地震參數以及井周圍地質構造。與地面地震相比，VSP資料的分辨率更高，可以得到更精確的時深轉換結果及速度模型，并為零相位子波分析提供支持。

前人針對VSP反射波資料使用基于射線理論的CDP疊加方法進行成像。Wyatt[1]首次提出VSP-CDP的概念，假設地層為常速水平層，當面對傾斜地層及復雜構造時，成像效果較差。Robert等[2]提出共中心深度點(CMD)疊加成像方法，基于水平反射界面假設，將共炮點道集數據轉換成CMD道集數據進行疊加成像，利用多道疊加提高了信噪比，但只適合反射界面傾角較小的情況。Smalley等[3]提出共側向點(CLP)疊加成像方法，假設在各向同性介質中共炮點道集滿足雙曲線時距方程，不適合橫向速度變化劇烈的地層。蔣振武等[4]研究了針對斜層的井間地震共深度點(DLCDP)成像算法。嚴又生等[5]闡述了非均勻介質井間地震VSP-CDP反射波疊加成像算法。石星等[6]應用共反射面元疊加成像方法對復雜構造成像，信噪比較常規疊加剖面更高。鄧金華等[7]提出了一種改進的共反射面元疊加方法，有效地避免了零井源距剖面中的同相軸相交現象。楊飛龍等[8]利用高斯束方法對VSP數據成像，提出了VSP高斯射線束法疊加成像算法?；诔Ｒ幧渚€追蹤理論的疊加成像方法主要解決VSP簡單波場成像問題，面對復雜構造時成像精度較低。

人們將基于偏移理論的成像方法引入VSP波場成像。Authur等[9]利用層析與相位屏法波動方程偏移聯合算法考察井中地震反射波成像。王華忠等[10]利用三維聲波方程實現了VSP 數據疊前深度偏移，并給出對應方程及其差分格式。宋煒等[11]使用Kirchhoff積分法進行井間地震反射波偏移成像。Liu等[12]以Kirchhoff偏移為基礎，使用波動路徑疊前偏移成像方法對井中地震數據成像。Fei 等[13]利用相移加插值(Phase-Shift Plus Interpolation，PSPI)法對二維、三維 VSP 數據偏移成像。劉詩竹[14]詳細研究了VSP數據的數值模擬和偏移成像方法。王維紅等[15]對Walkaway VSP逆時偏移進行了界面集成，推動了VSP 資料處理的實用化進程。蔡曉慧等[16]提出了基于自適應優化有限差分法的VSP 逆時偏移算法，在保證算法精度的前提下，提升了計算效率。楊繼東等[17]在高斯束偏移的基礎上，通過修改射線束傳播算子，提出了菲涅耳束偏移方法。楊飛龍[18]將高斯束偏移方法擴展到井中地震勘探，實現了井中地震高斯束疊前深度偏移成像方法。王沖等[19]利用Walkaway-VSP技術處理吐哈油田三塘湖盆地VSP資料，詳細刻畫了火山巖的形態和展布范圍。俞岱等[20]使用高斯束偏移方法試算井中地震模型及實際資料，通過對比成像結果與地面地震剖面驗證了方法的穩定性?；谄评碚摰某上窦夹g受觀測系統限制，導致成像剖面邊緣的“畫弧”現象嚴重。

在射線類偏移成像中，權函數對成像精度有一定影響。Liu等[21]研究了束偏移中的權函數，提出了適合Kirchhoff束偏移的疊加權函數。Sun等[22-23]在文獻[21]的基礎上改進權函數計算方法，提出了基于壓縮感知的Kirchhoff束偏移方法。若使用上述權函數對VSP數據進行VSP-CDP轉換，會出現波場能量失真，不能實現保幅疊加成像。

為此，本文使用動態射線追蹤方法模擬VSP復雜構造地震波場，通過交互對比正演波場與實際波場特征獲取準確的構造模型及速度場，最后使用正態分布權函數對VSP波場進行VSP-CDP轉換。前人在VSP-CDP轉換時僅根據地震波的傳播特征將共炮點道集數據中的每道數據歸位到相應的反射點處，本文使用正態分布權函數歸位VSP波場的同時，將深度—時間域的每一個采樣點轉換成反射點井源距—深度域的多個樣點，使反射點的覆蓋次數均勻，提高了VSP成像精度。

1 理論方法

1.1 動態射線追蹤

在直角坐標系下，波動方程可以表示為

(1)

式中：u為波場；v為傳播速度；t為傳播時間；x和z為笛卡爾坐標。

射線追蹤將波場分解為射線場，在已知射線參數的條件下，可選任一射線并建立相應的射線中心坐標系(s,n)(圖1)。

圖1 中心坐標系下射線傳播示意圖

對于任意一條射線，(s,n)是中心射線Ω上某一點的坐標，s代表地震波傳播的距離，n代表相鄰射線與Ω的垂直距離。坐標系的基矢量分別為與Ω相切的單位切向量t和與Ω垂直并指向Ω內側的單位法向量n

把式(1)變換到射線中心坐標系(s,n)，有

(2)

由式(2)可得到集中于射線中心鄰近的解。拋物型波動方程具有下列形式的時間調和解

(3)

式中ω為角頻率。式(3)表示沿中心射線進行積分。

對式(3)進行高頻近似并忽略ω的高階項可得u(s,n,t)的近似表達式

(4)

其中

(5)

Γ,s+vΓ2+v-2v,nn=0

(6)

(7)

式中：Γ=Γ(s)為位置的復值函數；A為相鄰射線與中心射線的權函數關系；Γ,s、A,s分別為Γ、A對s的一階偏導數；v,nn為v對n的二階偏導數。

式(6)為動態射線追蹤方程，式(7)為傳輸方程。由式(6)和式(7)解出Γ(s)和A(s)后可求出式(5)的W(s,v)，從而求得拋物型波動方程的解。

引入一新的復值函數q(s)，設

(8)

式中q,s為q(s)對s的一階偏導數。將式(8)代入式(6)，得到關于q的二階線性微分方程

vq,ss-v,sq,s+v,nnq=0

(9)

式中q,ss為q(s)對s的二階偏導數。再令q,s=vp，則式(9)可寫為

v(v,sq+vp,s)-v,s(vp)+v,nnq=0

(10)

式中v,s、p,s分別為v、p對s的一階偏導數。進一步化簡，得到一階微分方程組

(11)

將式(11)的第1式代入式(8)得

Γ=pq-1

(12)

將其變換到笛卡爾坐標系，得

(13)

(14)

式中：τ(s)為地震波旅行時；K(s)為地震波的波前曲率；Lh為地震波傳播到檢波點處的有效半寬度，與角頻率及動態射線參數有關。則式(13)變為

(15)

任意位置(s，t)的地震波場值為

(16)

式中：φ為射線入射的角度；uφ(s,n)為在(s,n)處的地震波場值。

1.2 正態分布權函數

常規的VSP-CDP疊加成像中，僅將地震反射波場轉換至相應地層的一個反射點處，難以滿足炮數少、構造復雜情形下的成像精度。本文從動態射線追蹤有效鄰域波場近似理論出發，在合成地震記錄時設檢波點處的波場為Lh范圍內所有射線的能量加權，因此進行VSP-CDP轉換時考慮有效范圍內所有射線能量的高斯加權，并根據權函數將檢波點的波場分解至Lh范圍內所有反射點處，大大提高了覆蓋次數。Liu等[21]認為： 射線束能量由中心射線往兩邊衰減，距中心射線越近，對應的振幅值越大；通過計算地震波走時發現，射線束內網格節點走時等相關信息是由泰勒近似求取的，距中心射線越近的點，相對誤差越小，精度越高，對應的權重越大。因此，Liu等[21]提出使用與距離相關的余弦平方窗函數作為成像公式中權函數主體，即

(17)

式中ns為相鄰射線與中心射線之間的距離。前人的研究結果以及動態射線有效鄰域波場近似理論表明，距檢波點越近的射線對能量的貢獻越大(圖2)。在VSP-CDP轉換時為了定量描述每條射線對檢波點能量的貢獻，假設檢波點處的波場值為1，則在Lh范圍內相鄰射線能量的加權等于1，并且距檢波點越近的射線對能量貢獻越大，符合正態概率密度函數分布。正態分布的概率密度函數與x軸圍成的面積為1，即將Lh范圍內的射線能量加權作為該檢波點的能量。根據正態概率密度函數的定義

(18)

式中μ、σ(σ>0)為兩個常數。σ越小，正態曲線越陡峭；σ越大，正態曲線越平坦。

參考正態曲線的形態及物理意義，可以將相鄰射線與中心射線權函數關系A近似表示為

(19)

當相鄰射線與中心射線重合時(ns=0)，由單一中心射線的波場能量表征該接收道的波場(常規射線追蹤方法)； 當nsLh時，相鄰射線的波場能量對該接收道的波場無貢獻。

圖2 中心射線與相鄰射線權函數關系示意圖

利用正態分布權函數將一道地震記錄進行VSP-CDP轉換至Lh范圍內有效反射點，那么有效反射點處能量的加權應為對應的檢波點波場能量。設某一道地震波場值為1，且該地震道的Lh為100m，利用權函數對該道地震波場進行VSP-CDP轉換，得到ns變化時反射點處波場特征示意圖(圖3)?？梢?，隨著相鄰射線與中心射線距離增大，地震波場值逐漸減小。若使用文獻[21]的權函數進行VSP-CDP轉換，則會出現有效射線的波場加權大于原始地震波場(圖3a)；本文提出的正態分布權函數使有效射線的波場加權等于原始地震波場(圖3b)，因此可保幅處理地震資料。

圖3 ns變化時反射點處波場特征示意圖

1.3 VSP-CDP轉換

1.3.1 常規VSP-CDP轉換原理

常速水平層狀介質VSP反射波時距曲線可以表示為(圖4)

(20)

式中：V為平均速度；x0為井源距；zR為檢波點深度；t為反射波旅行時。

VSP-CDP轉換的實質是把共炮點道集(或共檢波點道集)數據的每一個采樣點從深度—時間域變換到井源距—反射點時間(深度)域，形成共反射點道集數據。根據幾何關系(圖4)可知

(21)

式中(x1，z1)為反射點坐標。式(21)給出了常速水平地層情況下樣點從(z，t)域到(x，z)域的轉換公式。

圖5為單道數據的VSP-CDP轉換示意圖。由圖可見，根據上述VSP-CDP轉換方法可以計算到達檢波點R的兩個反射界面的反射點坐標(x1,z1)和(x2,z2)，分別對應T1和T2時刻，VSP-CDP轉換將共炮點道集記錄中T1、T2時刻之間的地震記錄分別歸位到反射點(x1,z1)、(x2,z2)處(若要轉換成井源距—反射點深度域，需要對T1時刻之前的地震記錄進行時深轉換，再進行反射點歸位)，這樣就形成了兩個反射界面情況下VSP單道數據的VSP-CDP轉換。同理，可完成單道數據所有反射界面及整個VSP多道數據的VSP-CDP轉換。

圖4 常速水平地層VSP-CDP轉換幾何示意圖

圖5 單道數據的VSP-CDP轉換示意圖

1.3.2 基于正態分布權函數的VSP-CDP轉換原理

在進行VSP-CDP轉換時，利用動態射線追蹤正演模擬計算反射點位置及檢波點有效鄰域范圍內射線對能量貢獻的正態分布權函數，并采用該權函數對波場分離后的VSP野外數據進行VSP-CDP轉換，將共炮點道集數據的一個樣點分解至共反射點道集數據的多個樣點，即將檢波點R的波場根據相鄰射線與中心射線正態分布權函數關系分解至Lh范圍內的所有射線對應的反射點處(圖6)。最后按照一定的間隔劃分面元，并對同一面元內的樣點進行疊加。對所有樣點進行上述操作，就得到最終的疊加成像剖面。與常規的VSP-CDP疊加成像方法相比，基于正態分布權函數的VSP-CDP方法有效增加了相同疊加面元內的反射點個數，進行VSP-CDP轉換后，有效擴大了反射點的照明范圍，大大改善了VSP勘探的橫向分辨率。

圖6 基于正態分布權函數的VSP-CDP轉換原理

2 模型試算

使用由表1參數表征的觀測系統對地質模型(圖7，模型參數如表2所示)進行VSP波場正演模擬，得到切除初至后的原始地震波場(圖8)。根據地質及測井資料建立初始速度場，并在同一觀測系統下進行VSP動態射線追蹤正演數值模擬(圖9)。調整速度及地質模型的過程是一個迭代過程，通過交互地質建模軟件微調地質模型，同時根據地震反射同相軸斜率修改地層速度，直至正演波場與原始波場吻合為止，圖10為第2炮數據原始波場與正演波場。使用基于正態分布權函數的VSP-CDP疊加成像方法進行反射P波和反射SV波成像，并選擇10m×10m疊加面元得到疊加剖面(圖11)?？梢?，基于正態分布權函數的VSP-CDP疊加成像方法對地下復雜構造及微小構造的成像精度較高。

表1 野外采集觀測系統參數

圖7 地質模型

表2 地質模型初始參數

圖8 切除初至后的原始地震波場

圖9 VSP動態射線追蹤正演波場

圖10 第2炮數據原始波場(黑色)與正演波場(紅色)

圖11 基于正態分布權函數的VSP-CDP疊加成像剖面

3 方法對比

為了驗證基于正態分布權函數的VSP-CDP疊加成像方法的效果，對地質模型(圖7)的局部區域使用常規射線追蹤的VSP-CDP疊加、基于正態分布權函數的VSP-CDP疊加方法成像，選擇10m×10m的橫向疊加面元得到疊加剖面(圖12)?？梢?，常規VSP-CDP疊加成像方法在靠近井旁構造區域及斷層面以下區域出現空道(圖12a)，基于正態分布權函數的VSP-CDP疊加成像方法沒有出現空道(圖12b)，表明前者的橫向分辨率低于后者。

圖12 不同成像方法效果對比

4 實際資料測試

為了驗證基于正態分布權函數的VSP-CDP疊加成像方法的效果，對銀額盆地哈日凹陷東部洼陷區Y5井非零井源距VSP資料進行成像。該VSP資料包括1炮零井源距數據、1炮非零井源距數據，圖13為Y5井非零井源距VSP三分量數據，圖14為Y5井非零井源距VSP三分量合成數據。根據研究區三維地面地震解釋結果建立實際地質模型(圖15)，通過零井源距VSP資料獲取井旁地層速度，然后使用VSP動態射線追蹤方法正演模擬，并將正演記錄與VSP上行反射縱波波場交互對比、修改速度模型，直至兩者吻合為止(圖16)。分別使用常規VSP-CDP疊加方法和基于正態分布權函數的VSP-CDP疊加方法對Y5井非零井源距VSP資料進行成像(圖17)?？梢?，無論是成像范圍還是成像精度，基于正態分布權函數的VSP-CDP疊加(圖17b)的成像結果都明顯優于常規VSP-CDP疊加(圖17a)。圖18為過井地面地震剖面與非零井源距VSP縱波時間域成像剖面。由圖可見：非零井源距VSP縱波時間域成像剖面(圖18b)與過井地面地震剖面(圖18a)波組特征一致，地層產狀相近；在2400～2700ms層段后者能量較弱，地層內幕結構信息模糊(圖18a)，前者的波組特征連續，成層性好，內幕信息豐富，分辨率高于后者(圖18b)。

圖13 Y5井非零井源距VSP三分量數據

圖14 Y5井非零井源距VSP三分量合成數據

圖15 實際地質模型

圖16 正演記錄與VSP上行反射縱波波場交互對比

圖17 非零井源距VSP縱波時間域成像剖面

圖18 過井地面地震剖面(a)與非零井源距VSP縱波時間域成像剖面(b)

5 結束語

基于動態射線追蹤有效鄰域波場近似理論，本文提出了一種非零井源距VSP資料疊加成像方法。通過研究正態分布的性質及特點，導出了基于正態分布疊加的權函數計算公式，并用于VSP疊加成像，改善了反射點覆蓋次數不均勻現象，提高了地震成像橫向分辨率。

模型試算及實際資料處理結果表明，基于正態分布權函數的VSP-CDP疊加成像方法能夠對復雜構造VSP資料精確成像，驗證了方法的有效性和穩定性。