近垂直定向分布裂縫的各向異性梯度地震反演方法研究

李春鵬,王 哲,李愛山,袁 峰

(1.中海石油 國際能源服務(北京)有限公司，北京 100028；2.中國石油大學(華東)，青島 266580；3.中海油能源發展股份有限公司 工程技術分公司，天津 300452)

0 前言

裂縫作為地殼中一種普遍的構造現象，廣泛存在于各類巖石中。到目前為止，已在砂巖、泥頁巖和碳酸鹽巖，甚至火成巖等各類巖石的裂縫性儲層中獲得了大量的工業油氣流。實際地層中裂縫呈一種至多種形態，包括近水平定向分布裂縫(VTI介質)、近垂直定向分布裂縫(HTI介質)、傾斜定向分布裂縫(TTI介質)和正交裂縫(OA介質)等。對于近垂直定向分布裂縫(HTI介質)，可以利用各向異性梯度進行預測，各向異性梯度越大，裂縫密度也越大[1-3]。

各向異性梯度可以通過方位地震數據反演得到[4]，寬方位比窄方位包含了更多的各向異性信息，實際生產中一般都利用寬方位地震數據反演各向異性參數[5-6]。不同于三維的疊后數據和四維的疊前數據，方位地震數據是五維數據，目前主要依據炮檢距向量片(OVT)技術處理寬方位地震數據，OVT技術是一種新穎的疊前數據的編排方式，基于OVT數據域處理可有效改善寬方位數據處理效果，且OVT域偏移結果含有豐富的方位各向異性信息，能實現寬方位數據的保真處理并提高疊前地震道集的質量，是寬方位三維地震數據的有效處理技術[7]。

地震反演的本質是求解線性方程組[8-9]，線性方程組能否得到良態的解，取決于它的系數矩陣的條件數。條件數越小，線性方程組的解越趨于良態；反之，越趨于病態。對于無約束地震反演方程，其系數矩陣是子波矩陣，然而子波矩陣的條件數特別大，所得到的反演結果是病態的。

為了解決這一問題，一般會對反演方程做加法，即加入高斯分布、柯西分布等數學約束，同時也會加入測井低頻約束[10-13]。事實證明，地震反演方程系數矩陣的條件數在加入數學約束和測井低頻約束之后顯著降低了，反演結果由病態變為良態。

實際工區一般都有測井縱、橫波數據，然而各向異性數據卻不常有。沒有各向異性數據也就沒有各向異性低頻約束，這樣就難以對各向異性反演方程做加法。各向異性梯度反演方程會因缺少測井各向異性低頻約束而難以降低反演方程條件數，從而難以得到準確的反演結果[9]。

各向異性梯度反演方程因缺少各向異性測井信息難以做加法，應對各向異性梯度反演方程進行其他改造。筆者基于這種思路，提出了對各向異性梯度反演方程做減法的方法，即穩定的各向異性梯度無井約束反演方法。該方法通過在反演方程中消除子波矩陣和降低反演方程的維數兩個減法操作，來降低反演方程系數矩陣的條件數，從而得到穩定的反演結果。

1 常規各向異性梯度反演方法

1.1 各向異性梯度的數學表征

高角度定向排列裂縫地層可以等效為HTI介質，Ruger[4]推導了上、下HTI介質對稱軸一致的方位縱波反射系數近似公式，Downton[14]在Ruger基礎上將該近似公式進行了線性化，即

(1)

圖1 地震反演流程圖Fig.1 Work flow of seismic inversion

(2)

通過式(1)可以建立反演方程，據此可以反演得到C和D，然后根據式(2)可以計算各向異性梯度。

1.2 常規各向異性梯度反演方程

方位地震數據的最小二乘反演方程為式(3)。

(WG)T(WG)m=(WG)TS

(3)

式中：S為方位地震數據;W表為子波矩陣;m=[ABCD]T為待反演參數;G為反射系數正演算子矩陣;R為反射系數。

(WG)T(WG)是反演方程的系數矩陣，該矩陣條件數極大，反演結果不穩定。一般都是對該系數矩陣做加法，通過加入高斯分布、柯西分布等數學約束和測井低頻約束來降低系數矩陣的條件數，以得到穩定的反演結果。但是實際工區有時會缺少各向異性梯度先驗分布信息和測井約束信息，因此難以通過常規反演方程得到真實的地下地層各向異性梯度。

2 各向異性梯度無約束反演方法

常規反演方程做加法的約束項的數學作用在于減小反演方程系數矩陣的條件數，以得到良態的反演結果。其實除了增加約束項，還有消除子波影響和降低反演方程維度的方法來降低反演方程系數矩陣的條件數，該方法分為兩步(圖1)。

1)改造反演方程，方法是子波移位到反演結果中。假設方位地震數據是方位反射系數，則方位地震數據正演方程如式(5)所示，式(3)和式(5)的最小二乘反演方程如式(6)和式(7)所示。

S=Gm'

(5)

m=[(WG)TWG]-1(WG)TS

(6)

m'=[GTG]-1GTS

(7)

2)疊后去子波，通過地震反褶積將第一步得到的wC和wD數據體轉換為C和D數據體，然后通過式(1)～式(5)可以計算各向異性梯度大小和各向異性梯度方向。但是常規反褶積方法是針對每一個采樣點計算的，計算量較大且不穩定。一般采用正則化算法改善反褶積算子，但是這種常規反褶積方法得到的反射系數不是稀疏的，不能真正地反映下地層界面的真實情況。筆者依據最少反射層地震反演方法[15]，先搜索目的層地下界面位置，只在搜索的界面位置處建立反褶積方程，這樣可以降低各向異性梯度地震記錄反褶積方程維數，進一步也會降低反褶積方程條件數，得到穩定的反演結果。

最少反射層稀疏脈沖反演的關鍵性技術是:首先利用模擬退火法尋找反射層位置，然后利用線性優化方法解目標函數，最終得到準確的反演結果[1,6]。

無約束各向異性梯度反演方法的關鍵是第一步子波移位，該步驟已經能得到待反演參數，即使反演參數含有子波項，也能體現地下地層的真實情況。另外該步驟簡單易操作，只需要把反演方程進行相應改造，即能得到穩定的反演結果。第二步是對第一步的補充，其關鍵是要能找到地下地層界面的位置，它受限于地震分辨率，一般適用于較厚的地層。實際操作中定義界面數目是同相軸的兩倍，能得到較好的結果。

綜上所述，無約束各向異性梯度反演方法和常規反演方法,對反演方程的改造存在較大的不同。無約束反演方法對反演方程做“減法”，減去子波項和非界面采樣點數，以降低反演方程條件數。而常規反演方法是對反演方程做“加法”，加入統計約束和測井低頻約束。雖然二者手段不同，但是目的都是增強反演方程的穩定性，以獲得準確的反演結果。但是這并不代表無約束反演方法可以代替常規反演方法，無約束反演得到的是相對結果，常規反演得到的是絕對結果且有更多的約束更符合地質規律。無約束反演更適用于缺少測井資料而難以提供低頻約束的情況，是地震反演的權宜之計，若資料允許，還是建議優先使用常規反演方法。

3 模型試算

為了驗證以上HTI介質各向異性梯度反演方法的有效性，引用四層三維模型[1-2,9](圖2(a))，其中第三層介質是近垂直定向分布的裂縫介質，第一層、第二層和第四層是各向同性介質。第三層裂縫介質的背景彈性參數和物性參數取自于A井的測井信息(圖3)，據此可以得到裂縫密度平面分布、裂縫走向平面分布和裂縫走向玫瑰圖(圖4)。根據以上模型數據可以建立HTI介質彈性矩陣，再根據HTI介質彈性波反射透射方程[16]計算方位反射系數[16]，模型正演采用相移30°后的50 Hz雷克子波，入射角分別是15°、25°和35°，方位角分別是0°、30°、60°、90°、120°、150°和180°。圖2(b)顯示了inline75，xline50處，入射角是35°，信噪比是10的方位縱波地震記錄，圖2(b)中最上面的同相軸表示模型第一層和第二層的各向同性界面的方位地震記錄，該地震記錄不隨方位角變化而變化；圖2(b)中最下面的兩個同相軸分別表示第三層高角度裂縫介質的頂和底的方位地震記錄，該地震記錄隨方位角變化而變化。

各向異性梯度“減法”反演可以分成兩步。

1)假設方位地震數據是方位反射系數，使得反演方程隱去子波項，則式(7)反演矩陣GTG條件數約為172。但若采用式(6)的反演方程，則反演矩陣(WG)T(WG)的條件數約為7.3×1017，遠大于式(7)計算的反演矩陣條件數，可見式(7)的反演方程穩定性更高。單道反演結果見圖5(a)，圖5中反演的自激自收反射地震記錄wA和各向異性梯度地震記錄wB、wC、wD與真實值吻合的較好，說明反演結果比較可信。另外，wC和wD平方根可以表示各向異性梯度地震記錄wBani大小，不難發現wC和wD平方根在各向同性界面處幾乎為零，在裂縫介質頂、底界面處有值，說明wBani能夠指示地層裂縫發育情況。

2)通過地震反褶積將第一步得到的wC和wD數據體轉換為C和D數據體，但是以上模型的反褶積算子wTw的條件數約為1.6×1018，反褶積方程極不穩定。因此筆者采用最少反射層假設反演C和D數據體，單道反演結果見圖5(b)，其中C和D的平方根可以表示各向異性梯度Bani大小，不難發現C和D的平方根在裂縫地層頂、底界面處有值，其余地層分界面處幾乎為零，說明最少反射層假設反演方法能準確體現裂縫地層的發育情況。

將穩定各向異性梯度反演方法推廣到整個四層介質模型中，比較圖6(a)和圖4(a)可以發現，反演的各向異性梯度大小和裂縫密度比較吻合，說明反演的各向異性梯度基本可以反映裂縫地層裂縫密度的平面展布；比較圖6(b)和圖4(b)可以發現,反演的各向異性梯度方向與裂縫地層裂縫走向基本一致或垂直，另外圖6(c)中各向異性梯度方向玫瑰圖中有兩個互相垂直的主方向，其中18°左右的主方向與圖4(c)中的裂縫走向主方向一致，因此各向異性梯度方向可以指示裂縫走向，但是存在90°不確定性，實際裂縫地層預測應用中,需要加入裂縫走向先驗信息以輔助識別地層真實裂縫走向。

圖2 四層介質模型Fig.2 Four layers model (a)示意圖；(b)信噪比10的方位地震記錄

圖3 L工區A井的測井曲線Fig.3 Well A logs of block L

圖4 第三層裂縫介質參數Fig.4 Fracture parameter of the third layer(a)裂縫密度；(b)裂縫走向;(c)走向玫瑰圖

圖5 反演結果(黑)與真實值(紅)對比Fig.5 Comparison between inversion result (red) and real result (black)(a)第一步；(b)第二步

圖6 反演結果Fig.6 Inversion result(a)各向異性梯度；(b)各向異性梯度方向；(c)方向玫瑰圖

圖7 不同方位疊后地震道集Fig.7 Post stack seismic of different azimuth(a)0°～30°;(b)30°～60°;(c)60°～90°;(d)90°～120°;(e)120°～150°;(f)150°～180°

圖8 裂縫預測結果Fig.8 Fracture prediction results(a)各向異性梯度；(b)各向異性梯度方向玫瑰圖；(c)測井結果

4 實際裂縫地層預測應用

利用某工區過A井的實際地震數據參與反演，該工區發育大段泥頁巖裂縫地層，工區內地震、測井資料比較完備，適合于利用方位道集進行裂縫型儲層預測研究。該方位地震數據的炮檢距和方位角分布比較均勻，屬于寬方位角地震數據，因此適合于利用方位道集進行裂縫型儲層預測研究。通過地震資料處理，得到了6個部分方位疊加道集，分別是0°～30°、30°～60°、60°～90°、90°～120°、120°～150°、150°～180°(圖7)，圖7中紅框內的裂縫發育區地震振幅隨著方位角變化而變化，說明方位地震對裂縫有較好的響應。

圖8是該工區S地層裂縫預測結果，已知A井在S地層附近存在裂縫，圖8中反演的S地層各向異性梯度在A井附近存在異常，這與已知認識吻合。另外，反演的各向異性梯度方向與測井的裂縫走向結果一致。以上表明，這里研究的各向異性梯度反演方法能得到準確的地層各向異性梯度。

筆者發展的各向異性梯度無約束反演方法，在第一步子波移位中就能得到較好的反演結果，該步簡單易操作且對地震資料分辨率要求不高。第二步去子波，利用模擬退火法搜索反射界面位置，對于調諧地層效果較差，因此適用于較厚的地層，實際應用需根據地震資料分辨率決定反演步驟。

5 結論

各向異性梯度對預測裂縫型儲層具有十分重要的意義，利用寬方位地震數據可以反演各向異性梯度。常規的地震反演是對反演方程做加法，但是對于沒有各向異性測井數據的工區，是難以增加測井低頻約束的。針對這種情況，可以利用筆者發展的消除子波和降低反演方程維度對反演方程做減法的方法，得到穩定的反演結果。該方法分為兩步：①子波移位簡單易操作能得到相對裂縫發育區;②去子波能得到絕對裂縫發育區，但操作較為復雜且對地震資料分辨率要求較高。模型試算和實際工區應用表明，該方法能夠得到準確的各向異性梯度，據此可以較為準確地預測高角度裂縫發育帶，進一步可以為鉆探裂縫型地層提供借鑒。