川西北雙魚石地區(qū)棲霞組地震資料優(yōu)化處理及裂縫預(yù)測技術(shù)應(yīng)用

于 豪，黃家強(qiáng)，蘭雪梅，劉軍迎

(1.中國石油勘探開發(fā)研究院,北京 100083；2.中國石油西南油氣田分公司勘探開發(fā)研究院,成都 610041)

裂縫是儲層的有效儲集空間,更是連通儲層和孔隙流體的運(yùn)移通道,直接影響油氣藏井網(wǎng)的部署。裂縫的建模、刻畫及展布規(guī)律是油氣藏開發(fā)的重要地質(zhì)依據(jù)[1]。

傳統(tǒng)的相干[2]、傾角[3]、曲率[4]等幾何結(jié)構(gòu)屬性是目前描述斷裂和裂縫發(fā)育帶的主要刻畫工具,這些屬性雖然反映了地層的不連續(xù)性,但同時也受到了地層巖性及其他非構(gòu)造因素的影響。雖然常規(guī)的地震隨機(jī)噪聲衰減技術(shù)可以在一定程度上消除這些影響,但是在壓制隨機(jī)噪聲的同時也損傷了有效信號,如斷點(diǎn)位置、地質(zhì)體邊緣等,不能有效識別地層邊界的真實(shí)形態(tài)[5]。

構(gòu)造導(dǎo)向?yàn)V波技術(shù)[6]具有壓制隨機(jī)噪聲、保護(hù)地層邊界信息的作用[7]?；谄⒎址匠痰母飨虍愋詳U(kuò)散方法[8]具有較高的穩(wěn)定性、計算效率以及較強(qiáng)的適應(yīng)性,被引入到對地震數(shù)據(jù)優(yōu)化處理的構(gòu)造導(dǎo)向?yàn)V波中。Perona等[9]在數(shù)字圖像增強(qiáng)技術(shù)領(lǐng)域提出了一種基于擴(kuò)散方程的噪聲衰減方法,即P-M模型,能夠保持邊緣有效信息和結(jié)構(gòu)特征。Weickert[10]通過引入結(jié)構(gòu)張量約束進(jìn)行多尺度非線性擴(kuò)散濾波,利用改進(jìn)的差分格式提高了計算效率。Fhemers等[11]在地震資料處理中提出了一種能夠約束擴(kuò)散濾波方向和強(qiáng)度、抑制地層邊界和斷點(diǎn)位置平滑的各向異性擴(kuò)散濾波技術(shù)。陳鳳等[12-13]提出了一種二維沿層濾波技術(shù)提高地震剖面的信噪比,并結(jié)合光流分析技術(shù)推廣到三維地震中。王旭松等[14]提出了一種一致增強(qiáng)性擴(kuò)散算法,最大限度地消除噪聲并保持邊界完整信息,并與孫夕平等[15]將非線性各向異性擴(kuò)散濾波技術(shù)應(yīng)用到二維地震剖面的去噪處理中。張爾華等[16]在三維地震資料中利用梯度結(jié)構(gòu)張量控制河道和斷層的邊緣特征,嚴(yán)哲等[17]將地震相干屬性和各向異性擴(kuò)散濾波相結(jié)合,以此作為斷層信息保護(hù)因子,隆雨辰等[18]將螞蟻?zhàn)粉櫯c曲率融合技術(shù)利用在川西地區(qū)斷層和裂縫刻畫中。近些年,地球物理學(xué)家們基于各向異性擴(kuò)散濾波研發(fā)了多種新方法,并取得了一定的應(yīng)用效果[19-24]。

另外,提高地震資料分辨率也能夠突出裂縫信息,特別是小尺度裂縫等細(xì)節(jié)特征,譜反演是提高地震分辨率的一種有效方法。Marfurt等[25]在薄窄河道的描述中提出了用滑動時窗的頻譜分析方法計算多種頻率相關(guān)屬性進(jìn)行儲層厚度定性解釋。Portniaguine等[26]認(rèn)為任何一個反射系數(shù)序列都可以分解為一個奇反射系數(shù)對和一個偶反射系數(shù)對,而地震響應(yīng)與這兩個反射系數(shù)對的大小有關(guān)。Puryear等[27]提出了Widess楔形模型理論,把反射系數(shù)序列分解成偶分量和奇分量,構(gòu)建了譜反演的目標(biāo)函數(shù)。Chopra等[28]由譜反演所得的反射系數(shù)算出了波阻抗剖面,并將其應(yīng)用到地層學(xué)解釋中。中國專家學(xué)者針對譜反演也開展了深入研究,例如基于隨機(jī)爬山反演方法[29]、基于模擬退火法[30]、基于最小二乘QR分解算法[31]、基于稀疏貝葉斯算法[32]、基于Cauchy條件約束算法[33]、基于Moore-Penrose算法[34]等。周家雄等[35]引入L-BFGS算法應(yīng)用于海上油田,劉萬金等[36]、司朝年等[37]和嚴(yán)皓等[38]分別將譜反演方法應(yīng)用于大慶油田、渭北油田和渤海油田,均得到了較好的應(yīng)用效果。

川西北雙魚石地區(qū)位于龍門山山前帶,受構(gòu)造活動影響,二疊系棲霞組發(fā)育北東-南西向斷裂,斷裂平面延伸范圍廣,強(qiáng)烈的構(gòu)造活動導(dǎo)致裂縫也比較發(fā)育。該地區(qū)儲層主要為一套碳酸鹽巖白云巖沉積序列,厚度較薄,差異較大,非均質(zhì)性較強(qiáng),多為孔隙-裂縫型儲層。目前該區(qū)已經(jīng)進(jìn)入開發(fā)階段,之前勘探階段的研究成果尚不能滿足開發(fā)的需要,急需厘清裂縫發(fā)育特征及展布規(guī)律。本文以雙魚石地區(qū)棲霞組為研究目標(biāo),首先對用于構(gòu)造解釋的疊后地震數(shù)據(jù)進(jìn)行了擴(kuò)散濾波和反射系數(shù)反演處理,在去除噪聲的同時,提高了資料的分辨率；通過多尺度屬性融合技術(shù)定性預(yù)測了裂縫發(fā)育特征；在只有疊后地震資料的前提下,利用神經(jīng)網(wǎng)絡(luò)和DFN(discrete fracture network)離散建模技術(shù)探索預(yù)測了裂縫的密度、長度和發(fā)育方向,為該區(qū)增儲上產(chǎn)提供了技術(shù)支撐。

1 研究區(qū)概況

1.1 地質(zhì)背景與地震資料概況

龍門山造山帶位于揚(yáng)子陸塊西緣,其逆沖變形始于晚三疊世末期的印支運(yùn)動,在整個燕山運(yùn)動減弱并持續(xù),在新生代再次強(qiáng)烈沖斷隆升,形成了龍門山褶皺沖斷帶-川西前陸盆地系統(tǒng)。龍門山?jīng)_斷帶自北西向南東發(fā)育4條大的主干斷裂,通常以北川-映秀斷裂為界將龍門山?jīng)_斷帶劃分為龍門山后山帶和前山帶。由于川西北地區(qū)構(gòu)造的形成演化直接受控于龍門山?jīng)_斷帶的發(fā)展,因此二者在現(xiàn)今構(gòu)造格局上顯示出很大程度的一致性[39]。雙魚石構(gòu)造帶位于川西斷褶帶西北緣,屬于中壩-雙魚石高帶的次一級正向構(gòu)造單元,北西面與射箭河-潼梓觀高帶相連,南東面與劍閣-梓潼坳陷帶相鄰。雙魚石地區(qū)構(gòu)造形態(tài)主要為北東向構(gòu)造,與龍門山近于平行,受龍門山推覆構(gòu)造控制,是印支、燕山、喜馬拉雅多期構(gòu)造運(yùn)動共同作用的結(jié)果(圖1)。

研究區(qū)為雙魚石三維區(qū),目的層為下二疊系棲霞組,主要為碳酸鹽巖臺地臺緣沉積,儲層主要發(fā)育為白云巖儲層。由于位于龍門山山前帶,所以構(gòu)造作用強(qiáng)烈,斷裂和裂縫較發(fā)育。研究發(fā)現(xiàn),雙魚石三維原始地震資料存在規(guī)則噪聲,從棲霞組頂界相干屬性切片上可以看出(圖2),地震資料存在與北東向斷裂斜交的規(guī)律性類指紋干擾,這會影響裂縫識別的精度,需要針對研究目標(biāo)開展解釋性處理,主要包括既能保持邊界特征又能去除規(guī)則噪聲的擴(kuò)散濾波處理,以及提高地震數(shù)據(jù)分辨能力的反射系數(shù)反演處理。

圖1 研究區(qū)地理位置及斷裂特征Fig.1 Location of the studied area and sketch map of the fault development

1.2 井中裂縫發(fā)育特征

根據(jù)研究區(qū)ST8井和ST10井在棲霞組的成像測井資料(圖3),得到井中裂縫發(fā)育特征。

ST8井裂縫密度為0.6,裂縫發(fā)育主方向?yàn)?0°。在7 332.5～7 347.5 m深度段,裂縫發(fā)育方向主要為0°～15°、60°～75°和150°～180°,傾角為40°～60°；在7 358.0～7 408.0 m深度段,裂縫發(fā)育方向主要為30°～45°,傾角在70°左右,裂縫類型為張開縫。

ST10井裂縫密度為1.2,裂縫發(fā)育主方向?yàn)?30°。在7 436.1～7 461.1 m深度段,裂縫發(fā)育方向主要為90°～165°,傾角在30°～50°；在7 461.1～7 486.1 m深度段,裂縫發(fā)育方向主要為150°～165°、195°～210°、330°～345°,傾角在30°～50°,裂縫類型為張開縫。

2 擴(kuò)散濾波

2.1 擴(kuò)散濾波原理

Perona和Malik[9]提出的基于偏微分方程的P-M擴(kuò)散濾波理論,實(shí)現(xiàn)了平行邊緣擴(kuò)散并保持邊緣特征,擴(kuò)散方程為

(1)

式(1)中：t為擴(kuò)散時間；div為散度算子；U為t時刻的擴(kuò)散濾波結(jié)果；U0表示t=0時的最初特征；D決定了擴(kuò)散方程的類型。當(dāng)D是標(biāo)量且為常數(shù)時,擴(kuò)散方程是線性各向同性方程；當(dāng)D是標(biāo)量且D=g(u)時,擴(kuò)散方程是非線性各向同性方程；當(dāng)D是張量且D=D(u)時,擴(kuò)散方程是非線性各向異性方程,可以沿著地層傾向較好的對邊緣進(jìn)行定位,保留不連續(xù)的構(gòu)造信息。

圖4 迭代次數(shù)對濾波效果的影響Fig.4 Influence of iteration umber on filtering effect

2.2 擴(kuò)散濾波效果對比

擴(kuò)散濾波關(guān)鍵參數(shù)主要有光滑步長、迭代次數(shù)和結(jié)構(gòu)張量維數(shù)三個參數(shù),其中迭代次數(shù)對處理結(jié)果影響較大。地震原始剖面在棲霞組主頻為30.5 Hz,帶寬為29 Hz,信噪比為0.97。對原始剖面進(jìn)行擴(kuò)散濾波,迭代5次后,棲霞組主頻為26.5 Hz,帶寬為29 Hz,信噪比為0.98；迭代10次后,棲霞組主頻為25.5 Hz,帶寬為29 Hz,信噪比為0.99(圖4)。對比分析發(fā)現(xiàn),濾波剖面相對原始剖面在去噪后變得相對干凈,迭代5次剖面在去噪的同時更能還原地質(zhì)特征,而迭代10剖面相對原始剖面在圖4中紅圈斷裂帶損失的有效信息較多,因此擴(kuò)散濾波最終選取迭代次數(shù)為5。從經(jīng)過擴(kuò)散濾波后的棲霞組頂界相干屬性切片上可以看到(圖5),地震資料消除了噪聲干擾,突出地震數(shù)據(jù)對斷層及裂縫的成像能力,斷裂條帶展布更加清晰。

3 反射系數(shù)反演

3.1 反射系數(shù)反演原理

反射系數(shù)反演是一種通過去除原始地震數(shù)據(jù)中的子波得到反射系數(shù)序列,從而在有限的帶寬內(nèi)加強(qiáng)有效信號的方法,能夠改善地震資料品質(zhì),提高分辨率。

在時間域內(nèi),假設(shè)一個兩層模型,當(dāng)分析點(diǎn)位于兩層之間時,地震反射系數(shù)表達(dá)式為

(2)

式(2)中：r1為第一層的反射系數(shù)；r2為第二層的反射系數(shù)；t為走時；T為兩層之間的時間厚度。對其進(jìn)行傅里葉變換可得:

g(f)=2recos(πfT)+i2rosin(πfT)

(3)

式(3)中：f為頻率；re和ro分別為反射系數(shù)對r1與r2的偶分量和奇分量。利用歐拉公式,令t=0,利用褶積計算可得兩層模型反射系數(shù)公式:

(4)

將其推廣到N層模型,并將其寫成矩陣形式,可得N層模型反射系數(shù)公式:

(5)

圖5 棲霞組頂界相干屬性(擴(kuò)散濾波后)Fig.5 The coherence of Qixia Formation top(after diffusion filtering)

圖6 不同主頻Ricker子波地震響應(yīng)特征Fig.6 Seismic response of different dominant frequency of Ricker wavelets

3.2 反射系數(shù)反演效果對比

在擴(kuò)散濾波的基礎(chǔ)上對地震資料進(jìn)行反射系數(shù)反演,為確保提頻的合理性,利用井資料進(jìn)行正演模擬,防止引入假頻或其他噪聲干擾?；趯?shí)際井資料,利用25～45 Hz的Ricker子波制作合成地震記錄(圖6),該井儲層在40 Hz主頻Ricker子波的時候,棲霞組頂界波峰下方波谷內(nèi)會出現(xiàn)一弱波峰反射,即該井儲層的地震響應(yīng),因此選用40 Hz作為提頻后的地震剖面主頻。

經(jīng)過反射系數(shù)反演前后地震剖面對比如圖7所示，可見提頻前剖面主頻32.5 Hz,提頻后剖面主頻39 Hz,主頻略有提高,頻帶在低頻端和高頻端都有明顯的拓寬,提頻后在棲霞組上部多出一個同相軸的地震響應(yīng),即薄儲層的響應(yīng)。從經(jīng)過反射系數(shù)反演后的棲霞組頂界曲率屬性切片上可以看到(圖8),提高分辨率以后地震資料反映的信息更加豐富,平面多出微小裂縫,除了可以反映區(qū)域分布的大斷裂和斷裂帶,還能表征沿斷裂分布的中小尺度裂縫。

圖7 反射系數(shù)反演前后地震剖面對比Fig.7 Seismic section before and after reflectivity inversion

圖8 棲霞組頂界相干屬性(反射系數(shù)反演后)Fig.8 The coherence of Qixia Formation top(after after reflectivity inversion)

4 裂縫預(yù)測技術(shù)應(yīng)用

4.1 多屬性融合裂縫表征

地震屬性是將地震數(shù)據(jù)經(jīng)過數(shù)學(xué)變換得出的有關(guān)地震波的幾何學(xué)、運(yùn)動學(xué)、動力學(xué)或統(tǒng)計學(xué)的特征反映,能夠反映地下豐富的地質(zhì)信息。在實(shí)際資料中,單一的地震屬性往往不能準(zhǔn)確全面地代表地下的真實(shí)面貌,需要多種屬性綜合分析,當(dāng)屬性之間存在一定程度的關(guān)聯(lián)時,可以進(jìn)行多屬性融合。

任意一種色彩都可以由紅、綠、藍(lán)各占一定的百分比形成,RGB模型就是一種由紅、綠、藍(lán)混合再生成其他色彩的模型,是一種在相同單位和值域范圍突出表現(xiàn)屬性的有效方法,因此在地震屬性分析中引入RGB色彩顯示技術(shù)進(jìn)行屬性融合。

針對不同尺度的裂縫,采用不同的幾何屬性如曲率、傾角和相干來表征。曲率是描述曲線上任一點(diǎn)的彎曲程度,在地震解釋中,曲率屬性是反映地震同向軸的彎曲程度,可以突出小尺度裂縫。受應(yīng)力場影響,研究區(qū)在地層彎曲的頂部更容易形成裂縫,利用最大正曲率表征該區(qū)裂縫發(fā)育情況。研究區(qū)棲霞組頂界最大正曲率范圍為-0.01～0.01(圖9)。

圖9 棲霞組頂界曲率屬性Fig.9 The curvature of Qixia Formation top

傾角是反映同向軸的傾向和走向,可以突出落差不到半個波場的隱蔽斷層,一般指示中等斷層及裂縫。研究區(qū)棲霞組頂界傾角范圍為0°～20°,其中10°～20°的傾角可以反映區(qū)域性斷裂及裂縫發(fā)育帶(圖10)。相干反映的是相似性問題,突出非一致性。地震資料的相鄰地震道會因?yàn)榇嬖跀鄬?、巖性突變、裂縫或特殊地質(zhì)體而引起波形變化,當(dāng)橫向變化大時,相干值??；當(dāng)橫向變化小時,相干值大。在傾角屬性的基礎(chǔ)上,采用帶地層傾角的第三代相干算法,計算全頻帶相干數(shù)據(jù)體,可以反映區(qū)域大斷裂和局部中小斷裂(圖8)。

圖10 棲霞組頂界傾角屬性Fig.10 The inclination of Qixia Formation top

圖11 棲霞組頂界RGB融合Fig.11 The RGB fusion of Qixia Formation top

將相干、傾角、曲率三種屬性進(jìn)行RGB融合(圖11),開展大、中、小不同尺度裂縫特征研究。通過多尺度屬性融合裂縫表征可見,雙魚石地區(qū)棲霞組頂界裂縫發(fā)育與區(qū)域應(yīng)力場及區(qū)域大斷裂密切相關(guān),裂縫多與斷裂伴生,沿北東-南西向呈條帶狀分布。由于工區(qū)西北部位于褶皺山前帶,擠壓應(yīng)力作用促使發(fā)育高密度不同尺度裂縫。

4.2 神經(jīng)網(wǎng)絡(luò)裂縫密度預(yù)測

BP網(wǎng)絡(luò)是一種具有三層或更多神經(jīng)元的神經(jīng)網(wǎng)絡(luò),包括輸入層、隱含層和輸出層。BP網(wǎng)絡(luò)的學(xué)習(xí)過程分為兩步:一是正向傳播,將輸入信息從輸入層傳到隱含層,經(jīng)過處理后將計算結(jié)果傳至輸出層；二是誤差反向傳播,如果經(jīng)正向傳播在輸出層沒有得到預(yù)想的結(jié)果,則把誤差信號沿原路徑返回,通過修改各層神經(jīng)元的連接權(quán)值并重新計算,減小誤差再輸出,迭代直至得到預(yù)期結(jié)果為止[40]。

神經(jīng)網(wǎng)絡(luò)反演方法不依賴于模型,具有較高橫向分辨率,適合強(qiáng)非均質(zhì)性裂縫預(yù)測。在已知井上裂縫發(fā)育密度的前提下,利用ST8井裂縫發(fā)育密度曲線和裂縫識別相關(guān)屬性如分頻相干、曲率等進(jìn)行神經(jīng)網(wǎng)絡(luò)學(xué)習(xí),得到裂縫發(fā)育密度體,用ST10井進(jìn)行檢驗(yàn)。預(yù)測結(jié)果表明(圖12),過ST8井剖面棲霞組裂縫發(fā)育稍低,密度為0.5～1；過ST10井剖面棲霞組裂縫發(fā)育較好,密度為1～1.5。預(yù)測結(jié)果與實(shí)際成像測井統(tǒng)計結(jié)果較吻合。

從棲霞組頂界神經(jīng)網(wǎng)絡(luò)裂縫密度預(yù)測平面圖中可以看出裂縫發(fā)育的密度分布范圍(圖13),研究區(qū)中部ST1井、ST8井周圍以及工區(qū)北部的山前帶均有裂縫發(fā)育富集區(qū),表明研究區(qū)裂縫發(fā)育受構(gòu)造活動影響較大。

圖12 神經(jīng)網(wǎng)絡(luò)裂縫密度預(yù)測剖面Fig.12 Seismic sections of fracture density by neural networks prediction

圖13 棲霞組頂界神經(jīng)網(wǎng)絡(luò)裂縫密度預(yù)測圖Fig.13 Fracture density of Qixia Formation top by Neural Networks prediction

4.3 DFN離散網(wǎng)格裂縫預(yù)測

通常情況下,測井反映的是厘米到米級的小尺度裂縫,地震反映的是百米到千米級的大尺度裂縫,而對于中尺度裂縫缺乏相應(yīng)的基礎(chǔ)資料。離散裂縫網(wǎng)格建模技術(shù)DFN可以有效解決中等尺度的裂縫刻畫問題。DFN是一種基于示性點(diǎn)過程的隨機(jī)建模方法,最早應(yīng)用于巖石工程等領(lǐng)域,后被用于油氣儲層裂縫建模。在建立油氣儲層DFN模型時,利用確定性裂縫和斷裂計算裂縫密度,點(diǎn)過程確定裂縫位置,示性過程確定點(diǎn)的屬性,如裂縫形狀、傾角、傾向、開度等屬性。

針對棲霞組頂界低頻相干屬性切片,通過基于Hessian矩陣的Sterger算法對裂縫進(jìn)行細(xì)化分割形成單像素圖像的二值圖,進(jìn)行裂縫識別圖像增強(qiáng)。通過基于Hough變換的裂縫矢量化技術(shù),提取確定的離散裂縫矢量,進(jìn)行定量化統(tǒng)計分析,提取發(fā)育長度和方向。同時結(jié)合神經(jīng)網(wǎng)絡(luò)裂縫預(yù)測提供的密度信息作為約束,進(jìn)行隨機(jī)離散建模,增強(qiáng)小尺度裂縫的表征,將疊后地震裂縫預(yù)測的尺度提高到單道道間距。從棲霞組頂界DFN離散建模裂縫預(yù)測圖可以看出(圖14),裂縫依然伴隨斷裂走向發(fā)育,北部山前帶發(fā)育密度較高。從裂縫地震預(yù)測及測井統(tǒng)計玫瑰圖(圖15)對比中可以看出,ST8井裂縫發(fā)育方向主要為北東向,ST10井裂縫發(fā)育方向主要為南東向,預(yù)測結(jié)果與成像測井統(tǒng)計結(jié)果吻合。

圖14 棲霞組頂界DFN離散網(wǎng)格裂縫預(yù)測Fig.14 Fracture prediction of Qixia Formation top by DFN

圖15 裂縫地震預(yù)測及測井統(tǒng)計玫瑰圖Fig.15 Results of seismic prediction and imaging log

5 結(jié)論

(1)川西北雙魚石地區(qū)棲霞組裂縫總體為北東-南西向伴隨斷裂走向呈條帶狀展布,局部發(fā)育為北西-南東向。北部山前帶由于構(gòu)造活動影響,裂縫發(fā)育密度高,為0.6以上；南部裂縫發(fā)育密度稍低,為0.4～0.5。

(2)各向異性擴(kuò)散濾波根據(jù)每次迭代出來的圖像梯度的大小來判斷圖像的邊緣,能較好地對邊緣進(jìn)行定位,且邊緣處的模糊程度較小。反射系數(shù)反演利用了薄層厚度與頻譜干涉的周期成倒數(shù)關(guān)系,將子波從原始地震數(shù)據(jù)中去除從而得到反射系數(shù)序列,在有限帶寬內(nèi)使有效信號得到加強(qiáng)。

(3)多屬性融合技術(shù)可以減少多解性,提高預(yù)測精度。RGB顯示技術(shù)是一種直觀的可視化工具,在地震裂縫屬性融合中,可以定性地展現(xiàn)更為豐富的地質(zhì)信息。在只有疊后地震數(shù)據(jù)的前提下,可以利用BP神經(jīng)網(wǎng)絡(luò)搭配DFN離散網(wǎng)格技術(shù)定量預(yù)測裂縫的密度、長度和發(fā)育方向,地震裂縫預(yù)測結(jié)果與成像測井資料相吻合,平面分布規(guī)律與斷裂展布方向較一致,說明方法具有可行性。