基于曲波域自適應導向濾波的海上采集腳印壓制方法

陳 平,秦 童,郭 軍,李少軒

(中海石油(中國)有限公司天津分公司,天津300459)

在海上和陸上地震勘探中廣泛存在采集腳印,給地震構造勘探和儲層地球物理研究帶來較大影響。海上拖纜的地震采集腳印與地震采集方向有高度的耦合性,其主要特點是沿著采集方向或某一個特定角度,淺層切片或方差數據體上出現周期性的振幅假象[1-2],淺層表現比深層表現更為明顯,有些會造成相位的移動,有些會造成時差的錯動,在剖面上可見類似垂直斷層的假象。采集腳印的存在對復雜斷裂活動區的斷層解釋帶來嚴重的干擾,當采集腳印影響范圍涉及到目的層時,會對砂體的尖滅點識別以及連通性判斷產生影響。采集腳印會直接影響勘探井位部署和開發方案設計,使海上油氣勘探的高效評價面臨較大挑戰。

采集腳印的壓制或消除方法分為3類[3]。第1類是優化采集設計的方法。王彥倉等[4]、碗學儉等[5]、駱宗強等[6]認為陸上采集腳印主要與觀測系統的滾動方式、炮線距和接收線距相關,減小滾動線數、增大橫縱比等可壓制采集腳印的影響;謝玉洪等[1]、王征等[2]、張旭東等[7]認為海上采集腳印還受非觀測系統影響,包括潮汐、水溫、含鹽度、水流方向等,采集時精確記錄上述參數有利于后續處理工作。第2類是改進疊前處理方法,如針對海上拖纜羽漂、轉圈造成采樣不規則的面元中心化技術,針對殘余時差的剩余靜校正技術[1],減小線距的OVT域內插法[8]和振幅均衡后的加權疊前偏移處理技術[9]。這兩類方法雖然能夠從照明、覆蓋次數和空間一致性的角度解決采集腳印問題,但耗時費力。第3類是進行疊后壓制的方法。此類方法較多,發展也較快,如CHOPRA等[10]提出的f-k域濾波法,GüLüNAY[11]提出的Kx-Ky濾波法,鄒振等[12]在其基礎上提出的拉普拉斯Kx-Ky濾波法,陳學華等[13]提出的一種基于自適應三維平穩小波變換的方法,蔡涵鵬等[14]提出的一種完備經驗模態分解的方法。此外,還發展出一些基于圖像處理的方法,如AL-BANNAGI等[15]提出的奇異值分解方法、魏天罡等[16]提出的Garbor濾波方法。這些疊后處理方法存在一定的適用性或局限性:f-k域濾波由于自身假頻和與有效信號疊置問題,或使斷層模糊、主要地震反射特征發生變化;同時由于海上拖纜周期性差,波數域濾波效果有時并不理想。而主成分分解或Garbor濾波會極大改變資料本身,只能用于解釋性處理;小波變換和經驗模態分解的基函數是各向同性的,分解的維度不足以緊致表征包含采集腳印及斷裂發育的復雜地震紋理結構,導致去除采集腳印后斷面、斷點等有效地震紋理受影響。

針對上述問題,本文提出基于曲波域自適應導向濾波的海上采集腳印壓制方法。利用曲波變換獨特的“似針狀”基函數和各向異性分解優勢,最優稀疏表征復雜的地震紋理,在海上地震資料全部時空范圍內,采集腳印被視為一種不規則噪聲,以地震優勢方向為導向,能在高維的曲波域與地震信號較好區分并得到壓制。

1 方法原理

小波變換類地震去噪方法是先對地震數據進行稀疏分解,即尋找到一個合適的稀疏變換基,在變換域中用較少的系數表示原始數據的主要特征,然后利用不同稀疏分解方法對細節刻畫的特點,設計相關閾值以消除特定類型的噪聲,達到提高信噪比的目的[17-18]。小波變換的基函數是各向同性的,以“點”為單位捕捉圖像特征,無法緊致表達存在復雜紋理結構的地震數據,而曲波變換具有很強的方向敏感性和獲取平滑輪廓的能力,在去噪效果和相對原始數據的保真度上有優勢[19-20]。曲波變換于20世紀90年代首次提出,第一代曲波變換又叫塊Ridgelet變換,但由于脊波基幾何形態不清晰,有時對具有奇異性的直線逼近效果不佳,隨后很快提出了第二代曲波基,其冗余度更小,且易于理解和實現,本文基于第二代曲波變換進行采集腳印的壓制。

信號在連續曲波域進行稀疏分解可寫成:

c(j,l,k)=〈f,φj,l,k〉

(1)

式中:c(j,l,k)為曲波系數;j,l,k分別為曲波函數的尺度、方向和位置參數;f為期望分解的信號;φj,l,k表示時空域的曲波函數。在頻率域實現曲波變換更為高效,用頻率域窗函數Uj代替φj,l,k,在頻率域極坐標(r,θ)下,定義頻率窗Uj為:

(2)

式中:[j/2]為j/2的整數部分;Uj為極坐標下的一個楔形區域,是一個個同心的圓環,符合各向異性的尺度特性,數值上正比于半徑窗W(r)和角窗V(θ)的乘積。依據這兩個窗函數可以對頻率域進行多尺度、多方向劃分,求同一尺度任意方向、位置上的φj,l,k可先將Uj(ω)反變換到時間域,然后對φj旋轉和平移。定義旋轉角度序列θl和平移參數k:

θl=2π·2[-j/2]·l,l=0,1,…,且0≤θl≤2π

(3a)

k=(k1,k2)∈Z2

(3b)

(4)

其中,Rθl表示旋轉θl弧度。則(1)式可寫為:

(5)

式中:f(x)為時空域信號。將時空域的內積轉換到頻率域,(5)式進一步寫成:

(6)

(7)

其中,Vj(Sθlω)是頻率域角窗,Sθl是一個剪切矩陣。離散曲波變換使用同心矩形環狀區域替代連續曲波變換中的環形區域,實現尺度和角度的劃分。在空間域,單個楔形窗對應的曲波基為“似針狀”的,在三維時空域垂直于脊的方向,單個楔形窗對應的曲波基是震蕩的,而平行于脊的方向是低通平滑的,因而能更好地描述地震紋理結構以及邊緣信息。

海上由非觀測系統導致的采集腳印與水流、風向、潮汐等相關,這些因素在全部采集時空范圍內可認為是隨機的,與拖纜方向一致的有限時空則為較強的相干噪聲,即采集腳印。假設垂直采集方向即聯絡線方向能看作連續時空的二維切片,采集腳印則近似等效為隨機噪聲,表現為時差錯動或相位變化,與垂直斷層尤為類似,一定程度上可認為是波形的“奇異點”。曲波變換的楔形基具有很強的非線性逼近能力,這些“奇異點”被分解到曲波域各個尺度、方向和位置上,與地震紋理的曲波域映射相比,其系數方向和能量更為發散,而地震有效信號則相對集中,因而在曲波域具備區分和壓制采集腳印的條件。

曲波按照分解的尺度可分為內層、中間層和最外層,越靠近內層頻率越低。內層用來表征信號的概貌,外層表征信號的細節和邊緣特性。除了前面描述的區別,采集腳印相對于地震的曲波系數還存在內層系數小、外層系數大,近垂直方向系數大、近水平方向系數小的特點。這要求壓制過程中曲波系數的閾值隨分解尺度(或層數)和方向同時變化,這樣才能在保留斷層及同相軸等有效信息的前提下壓制采集腳印。設計自適應閾值函數td為:

(8)

(9)

式中:cel為取整運算符,意為向上取整;len(C{s})為返回曲波系數C{s}在s尺度下分解方向個數;σ反映閾值整體的高低;α,β是經驗參數,反映尺度和方向對閾值影響的程度。(8)式中絕對值項用于返回邏輯值,控制優勢分解方向,默認為垂直和近垂直方向。由于曲波域1、3象限,2、4象限是中心對稱的,(9)式中用順時針180°區間變化的角度參數w0代替360°區間變化的角度參數w。

自適應閾值函數會隨分解尺度和分解方向變化而變化,即在同一尺度下,閾值會隨分解方向減小,而在同一方向上,閾值會隨尺度(或層數)增加而增加,同時增加優勢分解方向為約束條件,以滿足采集腳印壓制要求。在實際應用中,σ可先近似等于噪聲能級水平的1/10,固定σ依次測試參數α,β,為保證采集腳印的壓制效果,并有效保護較為復雜的斷層紋理,分解尺度s的指數權重α應小于2,分解方向w0的指數權重應不大于1。

將上面得到的閾值應用于歸一化后單位矩陣的曲波域系數E{s}{w0},有:

(10)

式中:L2為歐式范數運算符;nel為求矩陣元素個數運算符;C0{s}{w0}是單位矩陣的曲波域系數。采用硬閾值濾波,進而得到壓制后的曲波系數Ct{s}{w}:

Ct{s}{w}=C{s}{w}·(|C{s}{w}|≥td·
E{s}{w})

(11)

對壓制后的曲波系數進行反變換,便得到壓制采集腳印后的聯絡線地震剖面,對所有聯絡線重復上述操作,就能得到壓制采集腳印后的三維地震數據體。

2 模型試算

首先考慮單一水平界面的地震道模型,地震道數為500道,地震記錄長度為1000ms,采樣率為2ms,在500ms處設置反射界面,應用主頻為35Hz雷克子波進行褶積,并在第150道和第350道添加近似等效采集腳印的2～5個樣點時差。圖1為單一水平界面地震道模型去除采集腳印前、后的對比結果。

采用離散曲波變換進行分解,分解尺度為5。圖2 是去除采集腳印前、后地震道模型第1層至第4層曲波系數矩陣構成的圖像,逐層進行了歸一化處理。圖像內層表示大尺度分解的曲波系數,其中最內層是各向同性小波分解后的結果,圖像外層表示小尺度分解的曲波系數,各層都有4條邊,分別順時針存放3π/4至π/4,π/4至-π/4,-π/4至-3π/4,-3π/4至3π/4四個方向區間分解后的曲波系數。含采集腳印的地震道模型在非優勢方向出現許多微小的曲波系數,尺度越小或越靠近外層,曲波系數值越大,如圖2a 所示。直接濾除非優勢方向的全部曲波系數(圖2b),再反變換回時空域,采集腳印得到較為明顯的壓制,說明通過各向異性曲波稀疏分解,能夠較好分離和表征采集腳印與有效信號(圖1b)。

圖1 單一水平界面地震道模型去除采集腳印前(a)、后(b)的對比

圖2 單一水平界面地震道模型去除采集腳印前(a)、后(b)曲波系數對比

實際地層產狀是變化的,且構造更為復雜,設計接近實際地下反射結構的二維模型,如圖3a所示,水平方向為13000m,垂直方向為4000m,模型存在一條長期活動的邊界斷層及兩條伴生斷層,斷距隨斷層的規模和深度變化。用主頻30Hz子波正演得到二維模型的合成地震剖面,如圖3b所示。并給復雜模型的正演地震記錄添加部分采集腳印,如圖3c所示,可以看到明顯的同相軸抖動,容易解釋為“垂直斷層”。

圖3 復雜模型(a)、地震正演剖面(b)及含采集腳印地震剖面(c)

采用曲波變換對其進行5層分解,圖4是1～5層曲波域系數矩陣組成的圖像,可見主要地震能量集中于近垂直方向上,但其它分解方向除包含采集腳印外,也包含斷層等其它紋理信息(圖4a),因而不能直接濾除曲波系數。應用本文方法進行壓制后,有效信號的曲波域系數得到較好保留(圖4b),將壓制前、后的曲波系數矩陣直接相減,能看到系數差值在不同尺度和方向上是變化的(圖4c),且近垂直方向和與斷面有關的分解方向上系數差值較小,說明本文方法較為合理。

圖4 復雜模型正演地震記錄曲波域系數對比a 去除采集腳印前; b 去除采集腳印后; c 系數相減結果

圖5是壓制采集腳印前、后復雜模型地震剖面,壓制采集腳印后主要的地震反射特征基本不受影響,同相軸產狀變化自然,斷層及斷裂組合無明顯異常,采集腳印得到較好的壓制,資料品質得到改善。

圖5 壓制采集腳印前(a)、后(b)復雜模型的地震剖面

3 實際應用

曹妃甸某油田位于渤中凹陷西斜坡的構造脊上,斷層活動性強,有利于淺層油氣的運移和匯聚。油田及圍區發育北東東向和近東西向的斷層,斷層配置關系復雜,和淺層砂體互相耦合,形成一系列有利的構造-巖性圈閉。該區海上拖纜地震資料沿北偏西45°采集,拖纜作業時間近8個月,受冬季季風、漁業活動和潮汐等多種因素影響,加上處理時間較早,地震資料的采集腳印較為發育。圖6為曹妃甸某油田去除采集腳印前、后方差數據體切片,可見500ms、1000ms原始方差數據體切片均存在北—西向、與拖纜采集方向相關的采集腳印,1000ms方差數據體切片斷層“毛刺”現象較為明顯(圖6a、圖6c)。采集腳印的存在直接影響淺層斷裂系統的組合,也影響新近系淺水三角洲儲層的精細刻畫。

圖6 曹妃甸某油田去除采集腳印前、后的方差數據體切片a 500ms原始方差數據體切片; b 500ms去除采集腳印后方差數據體切片; c 1000ms原始方差數據體切片; d 1000ms去除采集腳印后方差數據體切片

采用本文方法進行采集腳印壓制,經過測試和實驗可知,當σ=0.004,α=1.6,β=0.5時,閾值函數能較好壓制本區采集腳印。如圖6b、圖6d所示,500ms方差數據體切片上消除了大部分采集腳印導致的異常相干能量,1000ms方差數據體切片上斷層“毛刺”現象得到一定減弱,斷裂特征更為明確,有利于厘清斷層的平面組合關系。

圖7是去除采集腳印前、后聯絡線地震剖面,可以看出,原始剖面上采集腳印從250ms延伸到950ms,淺層由于頻率較高,波組錯斷更為明顯,容易誤判為“垂直斷層”(圖7a),處理后在一定程度上消除了這種“垂直斷層”的解釋誤區,地震波組特征更為自然(圖7b),在差剖面中能看到壓制的采集腳印,且有效信號基本不受影響(圖7c)。

圖7 去除采集腳印前、后聯絡線地震剖面a 去除采集腳印前; b 去除采集腳印后; c 差剖面

基于新資料對全區已追蹤砂體重新排查后發現,部分砂體的尖滅位置和連通關系發生變化。圖8給出了去除采集腳印前、后潛力砂體的對比結果,可以看出,原始地震剖面上,975砂體高部位基本尖滅在采集腳印處的位置(圖8a),處理后的地震剖面上,975砂體往高部位還是連通的(圖8b),差剖面中對應位置的能量也較為集中在采集腳印附近(圖8c)。

圖8 去除采集腳印前、后潛力砂體的對比a 去除采集腳印前; b 去除采集腳印后; c 差剖面

綜合分析認為,該區采集腳印的存在不僅影響地震資料解釋結果,也影響構造-巖性控制下儲層的精細描述,對油田周邊挖潛有較大影響,采用有效的采集腳印壓制技術十分必要。

4 結論和建議

海上采集腳印由于周期性差和與斷面類似的紋理,常規方法無法有效壓制或改善原始地震資料品質。本文提出在垂直于采集方向可近似將采集腳印看作連續時空內的“奇異點”,利用曲波變換各向異性和地震紋理表達的優勢,將其映射到相比小波變換更高維度進行區分和壓制,并給出針對復雜地質構造的自適應導向濾波方法。模型和實際應用結果均證明了方法的有效性和適用性,能消除淺層大部分采集腳印造成的地質解釋假象。

本文方法除用于壓制采集腳印外,也可作為一種保邊緣提高信噪比的一般方法。實際應用中經驗參數對采集腳印壓制效果較為敏感,需按照一定順序不斷測試。另外,本文采用硬閾值的方法進行濾波,采用軟閾值濾波方法能進一步減少反變換后的“吉布斯”現象。