[發明專利]邊緣檢測低序級斷層識別方法

有效

• 申請號: 201810141769.3
• 申請日: 2018-02-11
• 公開（公告）號: CN108415077B
• 公開（公告）日: 2021-02-26
• 發明（設計）人: 蘇朝光;馬玉歌;韓宏偉;孫明江;孫金磊;揭景榮;劉啟亮;王釗
• 申請（專利權）人: 中國石油化工股份有限公司;中國石油化工股份有限公司勝利油田分公司物探研究院
• 主分類號: G01V1/30
• 分類號: G01V1/30
• 代理公司: 濟南日新專利代理事務所(普通合伙) 37224
• 代理人: 崔曉艷

1.邊緣檢測低序級斷層識別方法，其特征在于，該邊緣檢測低序級斷層識別方法包括：

步驟1，根據疊后地震數據體在時間空間域應用三維多級盲源分離與保構造濾波相結合的去噪技術壓制隨機噪聲，提高地震資料的信噪比；

步驟2，在時間域應用多道相干算法進行傾角搜索估計地層傾角值，在頻率域應用地層傾角在地震記錄上表現的時間延遲，對其進行傅里葉變換的時延特性來估計地層傾角；

步驟3，基于時間-頻率域混合傾角掃描，沿地層傾角在局部層拉平的基礎上計算分析窗內結構導向坎尼屬性；

步驟4，根據結構導向坎尼屬性數據體，通過連井剖面和沿層切片這些形式進行低序級斷層的地震識別；

在步驟1中，在時間空間域應用三維多級盲源分離與保構造濾波組合去噪技術壓制隨機噪聲；將獨立分量分析ICA技術與稀疏編碼Sparse Coding算法結合，構成三維多級盲源分離的ICA-SC算法；

在步驟1中，ICA-SC算法包括的步驟如下：

①采用ICA方法估計一個正交矩陣w，對地震數據x左乘w，使得ICA變換后的數據滿足稀疏分布；

②使用收縮函數通過極大似然估計的方法計算目標數據的各個稀疏分量

極大似然估計MLD表達式為：

M(x)＝sign(x)max(0,|x|-σ²|l'(x)|)

l(x)＝-ln p(x)

其中x為地震數據，s為尺度函數，記為s＝E{x²}，p(x)為拉普拉斯概率分布函數，l(x)為稀疏懲罰函數，σ²表示目標數據中噪音的方差；

由此可以得到拉普拉斯密度模型對應的收縮函數：

③對正交矩陣w求逆，得到

基于保構造平滑的基本理論，針對三維地震數據體中的樣點，選取一系列具有不同構造方位和形狀的分析窗，對分析窗中的數據采用三維高斯加權平均來取代常規的算數平均，進一步求取標準差并選取標準差最小的分析窗為最優時窗；

三維高斯加權平均求取均值公式如下：

式中x、y和t分別為橫測線、縱測線以及時間變量，u為分析窗口中的原始地震數據，σ_x，σ_y，σ_t是尺度參數，決定著平滑程度，其值越大，高斯濾波器的頻帶就越寬，平滑程度就越好；通過調節尺度參數，在邊緣特征模糊與平滑過程中由于噪聲所引起的過多的突變量之間取得折衷；

在步驟2中，地層傾角在地震記錄中表現為時間上的延遲，利用傅里葉變換的時延特性來計算地層傾角，假設在地面上接收到的某一道地震記錄為x(t),與它相鄰的地震記錄的時間延遲為t₀，則相鄰地震記錄道y(t)＝x(t-to),對x(t)的傅里葉變換由公式：

再依據傅里葉變換的時延特性公式：

求對應于不同角頻率的時間延遲公式：

上式中，X(ω)、Y(ω)為相鄰地震記錄的傅里葉變換，ω為角頻率，i為虛數單位；設有效頻段的起始值、終止值分別為ωstart，ωend，選擇地震記錄主頻附近的頻段，其中包含N個角頻率值，則這兩道地震記錄之間的時間延遲t₀公式：

即上式表明，最終得到的時間延遲值是各個不同頻率下對應的時延值的平均值。

2.根據權利要求1所述的邊緣檢測低序級斷層識別方法，其特征在于，在步驟3中，依據幾何關系進行傾角掃描，計算出傾角、方位角兩個數據體，然后轉換成沿著X，Y方向的延遲時間，生成結構導向屬性體。

3.根據權利要求1所述的邊緣檢測低序級斷層識別方法，其特征在于，在步驟3中，把數據體內的每一點都當成是目標點，分別提取一個大小為(2n+1)×(2n+1)的局部窗口；然后沿著地層法線對稱地提取K個局部窗口，窗口總數為2K+1個。

4.根據權利要求3所述的邊緣檢測低序級斷層識別方法，其特征在于，在步驟3中，利用坎尼算子模板分別求取每個局部窗口的不同方向的坎尼屬性，并對不同方向的屬性值進行平方求和，然后沿著法線方向對所有局部窗口求得的坎尼屬性作和，最后按能量對結果進行歸一化處理。

5.根據權利要求4所述的邊緣檢測低序級斷層識別方法，其特征在于，在步驟3中，采用以下步驟計算結構導向坎尼屬性：

(a)設二維高斯濾波函數為：

G(x，y)＝(1/2πσ²)exp(-(x²+y²)/2σ²)

式中：x，y表示該點在坐標軸上的坐標，σ為高斯標準差，它決定了平滑程度，把梯度矢量▽G的2個濾波卷積模板分解為2個一維的行列濾波器：

式中：kx，ky表示x，y方向上的空間波數；

這兩個卷積模板分別與原始地震數據f(x,y)運算：

即：

f_s＝f(x，y)×G(x，y)

(b)經過步驟(a)處理的像素點I(x，y)的方向和強度是基于2×2模板的一次微分的有限差分得到的，x、y的偏微分P_x[i，j]、P_y[i，j]為：

P_x[i，j]＝(I[i，j+1]-I[i，j]+I[i+1，j+1]-I[i+1，j])/2

P_y[i，j]＝(I[i，j]-I[i+1，j]+I[i，j+1]-I[i+1，j+1])/2

其中：I(x，y)為像素點坐標(x、y)上的灰度值；

利用坐標轉換可以把梯度的方向和強度變換到極坐標下，以二階范數為基礎來表示如下：

梯度的方向表示成：

θ[i，j]＝arctan(P_y[i，j]/P_x[i，j])

(c)利用非極大抑制處理經過步驟二處理的圖像，設像素點為I(i，j)，而其梯度方向設為0(i，j)，梯度方向0(i，j)為四個區域中方向的一個，將像素點I(i，j)分別與四個不同方向的相鄰數據點比較大小從而計算局部區域內的極大值；若I(i，j)的梯度方向在方向4區中，則分別將其梯度強度M(i，j)和4區中右上角、左下角的兩個點的梯度強度比較大小，如果不是最大值則將I(i，j)的M(i，j)置為0，相反地就將其置為1；

(d)分別設置兩個閥值V1和V2給非極大抑制后的圖像數據，像素點梯度幅值小于V1的將其值置成0，從而生成圖像A；接著再運算尋找到幅值小于V2的置成0，生成圖像B；由于V2V1，圖像B干擾噪音相對很小，可是同時也導致了邊緣丟失；而相應的圖像A邊緣信息比較豐富；把圖像A的邊緣信息填補到圖像B中以連接邊緣線，這就是雙閥值的確定與邊緣鏈接的過程。