楊 勝，滕奇志

（四川大學(xué) 電子信息學(xué)院圖像信息研究所，成都610065）

0 引 言

隨著地質(zhì)石油勘探開發(fā)的逐漸深入，石油與天然氣勘探方向開始逐漸由淺轉(zhuǎn)向深部、由常規(guī)油氣藏轉(zhuǎn)向特殊油氣藏。特別是裂縫性油氣藏，成為地質(zhì)石油勘探方向的主要研究目標(biāo)［1］。在裂縫性油氣藏中，研究裂縫的發(fā)育程度和分布范圍，對(duì)石油和天然氣開采具有重要的意義。通過圖像處理技術(shù)，對(duì)巖心CT圖像進(jìn)行處理和分析，是研究巖心裂縫的一種較好的方法。為了計(jì)算每條裂縫的相關(guān)參數(shù)，首先需要將裂縫進(jìn)行提取，再進(jìn)行統(tǒng)計(jì)分析相關(guān)參數(shù)。但有些裂縫存在交叉的情況，需要通過算法自動(dòng)地將交叉裂縫分離。圖1展示了裂縫CT圖像經(jīng)過閾值分割后，存在交叉裂縫的情況。

從圖1（c）看來，單個(gè)交叉裂縫二值圖像在局部上具有直線特征，可以通過提取交叉裂縫的直線部分來將交叉裂縫進(jìn)行分離。所以針對(duì)交叉裂縫分離的問題可以轉(zhuǎn)換成圖像直線提取的問題。20世紀(jì)60年代，F(xiàn)reeman提出了直線檢測(cè)準(zhǔn)則［2］。通過這個(gè)準(zhǔn)則，可以檢測(cè)出圖像中的直線部分。除此之外，有些學(xué)者在此基礎(chǔ)上進(jìn)一步引入了直線檢測(cè)方向優(yōu)先級(jí)的概念［3-4］。方向優(yōu)先級(jí)的引進(jìn)，提高了直線檢測(cè)過程中平滑程度，可以更有利于處理直線交叉的情況。本文基于Freeman鏈碼優(yōu)先級(jí)直線提取算法的思想，對(duì)巖心二維CT序列圖像的交叉裂縫進(jìn)行分離，并在算法基礎(chǔ)上引進(jìn)了概率優(yōu)先準(zhǔn)則，使算法更好地適用于交叉裂縫的情況。

圖1 二值化處理裂縫CT圖像Fig.1 Crack CT image for binarization

1 Freeman鏈碼交叉裂縫分割算法描述

Freeman鏈碼直線提取算法，主要是檢測(cè)圖像的直線部分。首先需要對(duì)裂縫二值圖進(jìn)行細(xì)化處理，得到裂縫骨架［5-6］。理想情況下，經(jīng)過圖像細(xì)化后得到的骨架是單像素寬度的，單像素寬的裂縫骨架更有利于裂縫骨架的特征提?。?-8］。但在實(shí)際情況中，在骨架的轉(zhuǎn)折處或者交點(diǎn)處，并不能完全實(shí)現(xiàn)單像素。所以，要先對(duì)骨架圖像進(jìn)行單像素處理［9］，然后再對(duì)裂縫骨架進(jìn)行Freeman鏈碼直線提取。由于巖心裂縫只在局部上具有直線特征，因此在對(duì)裂縫骨架進(jìn)行Freeman鏈碼直線提取后，會(huì)造成裂縫骨架的斷開。

本文在提取裂縫骨架直線的基礎(chǔ)上，考慮了裂縫骨架端點(diǎn)附近的梯度信息。根據(jù)骨架端點(diǎn)附近的梯度信息，可以使得斷裂的骨架直線合并。合并后便可以得到一條完整的裂縫骨架直線，由此便可以將交叉裂縫中單個(gè)裂縫的骨架分別提取出來。得到裂縫的骨架數(shù)據(jù)后，可以通過骨架平移操作，將完整的裂縫目標(biāo)提取出來。具體算法流程如圖2所示。

圖2 算法流程Fig.2 Algorithm flow

1.1 Freeman鏈碼直線提取

Freeman鏈碼是用線段的起點(diǎn)加上由八個(gè)方向符所構(gòu)成的一組數(shù)列，數(shù)列中的值代表著八鄰域中兩像素的相對(duì)位置。對(duì)于每個(gè)像素點(diǎn)來說，其相鄰的點(diǎn)都可以用八方向符中某一個(gè)值表示。八方向鏈碼如圖3所示。由此可知，裂縫骨架圖像中，任何一組骨架曲線都可以由一組鏈碼串表示。Freeman直線提取算法，從裂縫骨架的一個(gè)點(diǎn)跟蹤到另一個(gè)點(diǎn)，一定會(huì)對(duì)應(yīng)一組鏈碼序列。

圖3 八鄰域鏈碼Fig.3 Chain code of eight neighborhood

Freeman根據(jù)直線在離散平面的特性，提出了直線在離散平面時(shí)鏈碼應(yīng)該遵循的三條準(zhǔn)則［2］（簡(jiǎn)稱Freeman準(zhǔn)則）：

（1）鏈碼中至多只有2個(gè)方向碼，其中一個(gè)為主方向，且其值相差1。

（2）若有兩個(gè)方向碼，則次級(jí)方向碼必單個(gè)出現(xiàn)，不能連續(xù)出現(xiàn)。

（3）次級(jí)方向碼總是盡可能均勻地出現(xiàn)在鏈碼中。

由于裂縫骨架在局部上才完全符合直線特征，故完全依照Freeman準(zhǔn)則來提取裂縫是不可行的。本文算法根據(jù)裂縫骨架特征，對(duì)Freeman準(zhǔn)則做出的部分改進(jìn)如下：

裂縫在全局上的趨勢(shì)走向并不像純粹的直線那樣單調(diào)，2個(gè)方向碼并不能很好地表征一個(gè)裂縫在全局上的趨勢(shì)走向，故本文算法在Freeman準(zhǔn)則基礎(chǔ)上的至多兩個(gè)方向碼，額外增加了一個(gè)方向碼。其中一個(gè)為主方向且另兩個(gè)方向碼和主方向的值相差1。另外，在此基礎(chǔ)上增加了每個(gè)方向碼出現(xiàn)的概率信息，在對(duì)裂縫骨架每個(gè)點(diǎn)進(jìn)行迭代時(shí)，更新目前鏈碼中所有方向碼的概率。

統(tǒng)計(jì)鏈碼中每個(gè)方向碼的概率的好處是：

（1）由于裂縫趨勢(shì)走向的不確定性，鏈碼的主方向在每次迭代過程都可能會(huì)改變，利用概率信息可以快速得到當(dāng)前迭代時(shí)的主方向。

（2）當(dāng)鏈碼有3個(gè)方向碼時(shí)，在某一輪的迭代跟蹤過程中，發(fā)現(xiàn)3個(gè)方向碼都存在目標(biāo)點(diǎn)，這時(shí)選擇繼續(xù)跟蹤概率最大的那個(gè)方向碼，保證跟蹤過程盡可能平滑。

如對(duì)某個(gè)直線L進(jìn)行跟蹤。當(dāng)搜索到直線L的第一個(gè)像素點(diǎn)時(shí)，按照以下方式進(jìn)行迭代：

（1）當(dāng)前鏈碼為空時(shí)，則按照0～7的方向，逆時(shí)針?biāo)阉鞯谝粋€(gè)鏈碼值，如圖4（a）所示。

（2）當(dāng)前鏈碼中只有一個(gè)方向D1時(shí)，先判斷D1方向是否可以繼續(xù)跟蹤。若不能，則尋找另一個(gè)鏈碼方向D2，如圖4（b）所示。

（3）當(dāng)前的鏈碼中有兩個(gè)方向D1和D2時(shí)，首先跟蹤概率較大的方向，若失敗則繼續(xù)跟蹤另一個(gè)方向，若兩個(gè)方向全部跟蹤失敗，尋找新的方向D3，如圖4（c）所示。

（4）當(dāng)前的鏈碼中有3個(gè)方向D1、D2和D3時(shí)，按照從大到小的概率依次跟蹤，如圖4（d）所示。

圖中虛線的含義是，當(dāng)實(shí)線的方向全部不能繼續(xù)跟蹤時(shí)可以選擇的新方向。

圖4 跟蹤示意圖Fig.4 Tracing diagram

若圖像中分布有很多噪聲，F(xiàn)reeman直線檢測(cè)算法在跟蹤的過程中，也可以直接去除圖像中的毛刺和噪聲［10］。通過設(shè)置一個(gè)鏈碼閾值長(zhǎng)度，一輪迭代完成后，判斷當(dāng)前鏈碼串的長(zhǎng)度是否大于該閾值長(zhǎng)度，大于閾值長(zhǎng)度則說明該鏈碼串是合格的鏈碼串，而小于閾值長(zhǎng)度，則說明這個(gè)鏈碼串是噪聲應(yīng)當(dāng)舍棄，繼續(xù)掃描下一個(gè)骨架點(diǎn)，開始下一個(gè)鏈碼串的跟蹤過程。每次迭代完成后，都記錄已跟蹤過的點(diǎn)，若重新跟蹤到標(biāo)記跟蹤過的點(diǎn)則直接跳過，直至所有的骨架點(diǎn)都掃描完畢，結(jié)束本次跟蹤過程。實(shí)驗(yàn)結(jié)果如圖5所示。

圖5 Freeman鏈碼直線提取Fig.5 Line extraction for Freeman chain code

1.2 裂縫骨架合并

由上述可知，在對(duì)裂縫骨架進(jìn)行鏈碼跟蹤之后，會(huì)出現(xiàn)裂縫骨架斷裂的情況，這時(shí)需進(jìn)行骨架合并計(jì)算。裂縫骨架是否可以合并和裂縫骨架端點(diǎn)之間的距離及骨架的偏轉(zhuǎn)角有關(guān)。設(shè)第L k條骨架的首尾端點(diǎn)為L(zhǎng) k（x1，y1），L k（xm，y m），根據(jù)公式（1）：

可獲得裂縫骨架的偏轉(zhuǎn)角θL k。判斷骨架L k和L k+i是否可以合并的充要條件是：偏轉(zhuǎn)角滿足公式（2），骨架端點(diǎn)間的距離滿足公式（3）：

其中，di s t函數(shù)代表兩坐標(biāo)點(diǎn)間的距離，α、β分別代表骨架偏轉(zhuǎn)角閾值和骨架端點(diǎn)距離閾值。當(dāng)裂縫骨架經(jīng)過跟蹤算法處理后的兩條裂縫骨架滿足裂縫合并的充要條件時(shí)，便可將二者進(jìn)行合并，具體算法如下：

（1）搜索到還未處理過的裂縫骨架L k，轉(zhuǎn)步驟（2）；若檢測(cè)到全部裂縫骨架已經(jīng)處理，則轉(zhuǎn)到步驟（3）。

（2）判斷裂縫骨架L k和L k+i是否滿足合并的充要條件。若滿足，則將兩條裂縫骨架合并為一條裂縫骨架L k，并將L k+i設(shè)置為另一條未處理過的裂縫骨架，重復(fù)步驟（2），否則轉(zhuǎn)到步驟（1）。

（3）結(jié)束裂縫骨架合并算法。

在上節(jié)裂縫骨架完成Freeman鏈碼直線提取的基礎(chǔ)上，進(jìn)行骨架合并。實(shí)驗(yàn)結(jié)果如圖6所示。

圖6 骨架合并Fig.6 Skeleton merging

1.3 裂縫完整目標(biāo)提取

由圖6結(jié)果可知，最終得到是一條線性單像素寬度的裂縫骨架，并不能完全覆蓋真實(shí)的裂縫區(qū)域。所以，需要通過一定的算法提取整個(gè)裂縫目標(biāo)。通過對(duì)裂縫骨架平移后，與裂縫二值圖像求交集操作，可以將完整的裂縫目標(biāo)提取出來。裂縫骨架平移方向上的選擇和裂縫骨架鏈碼的主方向呈90°關(guān)系，并且由于細(xì)化的裂縫骨架處于裂縫二值圖的中軸，所以骨架的平移方向是兩個(gè)互補(bǔ)的方向。算法實(shí)現(xiàn)過程如下：

（1）將方向平移間隔變量n1和n2初始化為1。

（2）掃描一條裂縫骨架L k，將L k的坐標(biāo)數(shù)據(jù)拷貝到P。根據(jù)裂縫骨架鏈碼的主方向求出平移方向，轉(zhuǎn)到步驟（3），若完成所有裂縫骨架平移后，轉(zhuǎn)到步驟（5）。

（3）將P所有坐標(biāo)拷貝后，向平移方向平移n1個(gè)間隔，并且與裂縫細(xì)化之前的二值圖相比較。若重合度達(dá)80%，則將平移后的坐標(biāo)數(shù)據(jù)更新到裂縫骨架L k中，n1自加1，重復(fù)步驟（3），否則轉(zhuǎn)到步驟（4）。

（4）將P中所有坐標(biāo)拷貝后，向另一方向平移n2個(gè)間隔，并與裂縫二值圖相比較。若重合度達(dá)80%，將平移后的坐標(biāo)數(shù)據(jù)更新到裂縫骨架L k，n2自加1，重復(fù)步驟（4），否則轉(zhuǎn)到步驟（1）。

（5）結(jié)束裂縫骨架平移算法。

通過上節(jié)的裂縫骨架提取裂縫，完整目標(biāo)結(jié)果如圖7所示。

圖7 裂縫完整目標(biāo)提取Fig.7 Extraction of crack complete target

2 實(shí)驗(yàn)結(jié)果

為了驗(yàn)證本文算法的正確性，取4組真實(shí)的巖心CT圖像作為實(shí)驗(yàn)數(shù)據(jù)，實(shí)驗(yàn)結(jié)果如圖8所示。其中圖8（a）表示在真實(shí)環(huán)境下存在交叉裂縫的巖心CT圖像；圖8（b）是將圖8（a）經(jīng)過閾值分割、去噪等形態(tài)學(xué)操作提取裂縫的結(jié)果；圖8（c）是由本文算法對(duì)圖8（b）處理的結(jié)果。

圖8 實(shí)驗(yàn)結(jié)果Fig.8 Experimental result

由上述實(shí)驗(yàn)結(jié)果可知，本文算法能夠在不影響裂縫本身特征的情況下將交叉裂縫進(jìn)行分離，并且分離結(jié)果符合現(xiàn)實(shí)預(yù)期結(jié)果，為后續(xù)裂縫統(tǒng)計(jì)參數(shù)分析提供更可靠的數(shù)據(jù)支撐。

3 結(jié)束語(yǔ)

傳統(tǒng)裂縫提取算法在遇到交叉裂縫的情況下，無法進(jìn)行有效的處理。為此本文提出了一種基于Freeman鏈碼的交叉裂縫分離算法。該算法加入了概率統(tǒng)計(jì)的因素，使追蹤裂縫的過程更加平滑；使用距離約束達(dá)到去除圖像噪聲的效果；對(duì)分割錯(cuò)誤的裂縫進(jìn)行了合并處理，減少了過分割情況的出現(xiàn)。實(shí)驗(yàn)結(jié)果表明，基于Freeman鏈碼的交叉裂縫算法對(duì)交叉裂縫的分割有著良好的效果，且算法計(jì)算量小、抗噪性強(qiáng)，適于圖像的實(shí)時(shí)處理，并為巖心地質(zhì)裂縫的參數(shù)計(jì)算提供了更加精密的計(jì)算方式。