段友祥,李寧寧,孫歧峰,李 鈺
(中國石油大學(xué)(華東) 計算機(jī)與通信工程學(xué)院,山東 青島 266580)
地質(zhì)剖面圖是用規(guī)定的符號、花紋和顏色,按一定的比例,沿一定的方向,表示一定的距離內(nèi)地下一定深度內(nèi)地質(zhì)現(xiàn)象的圖件,一般要表示出地層分層、斷層、巖性等地層結(jié)構(gòu)和巖體屬性特征。它形象直觀地表達(dá)了地層的結(jié)構(gòu)構(gòu)造和地層的沉積規(guī)律,是地層在垂向上最直觀、有效的表達(dá)方式。油藏地質(zhì)剖面圖能為系統(tǒng)分析區(qū)域或局部的油藏地質(zhì)條件、正確指導(dǎo)油藏資源開發(fā)利用和優(yōu)化管理提供依據(jù)。
現(xiàn)在很多老油田已經(jīng)進(jìn)入了后期開發(fā),小油藏、隱蔽油藏等復(fù)雜油藏成為主要的開發(fā)對象,運用的技術(shù)手段也越來越多樣化,要求的基礎(chǔ)地質(zhì)數(shù)據(jù)資料也越來越多。但是,由于受早期條件的限制,這些老油田的很多地質(zhì)資料不全或不符合現(xiàn)在開發(fā)技術(shù)手段的要求。例如,作為重要地質(zhì)成果圖件的地質(zhì)剖面圖,過去大多采用手工繪制(很容易轉(zhuǎn)化為位圖圖像)或軟件繪制位圖圖像(又稱為柵格圖像)。它們雖然能夠逼真展現(xiàn)地層剖面,但是難以運用現(xiàn)代計算機(jī)技術(shù)進(jìn)行處理,圖像的數(shù)據(jù)量大,不能進(jìn)行保真縮放,難以滿足現(xiàn)代的鉆井導(dǎo)向、油藏精確描述識別等先進(jìn)技術(shù)的要求。因此,對過去陳舊的地質(zhì)資料進(jìn)行深度加工和處理是老油田有效開發(fā)的重要基礎(chǔ),矢量化就是其中重要的內(nèi)容之一。矢量化是將油藏地質(zhì)剖面圖等柵格圖像轉(zhuǎn)換為矢量圖像的過程,其實質(zhì)是將圖像數(shù)據(jù)轉(zhuǎn)換為矢量圖形數(shù)據(jù),從而更方便地對圖形進(jìn)行處理,并精確地表達(dá)空間對象的位置、長度等信息[1]。
隨著計算機(jī)技術(shù)的發(fā)展以及圖像處理技術(shù)應(yīng)用的深入,圖像矢量化對象由早期的二值位圖(即黑白圖片)變?yōu)楝F(xiàn)在的彩色圖像,矢量化算法也在不斷發(fā)展。圖像矢量化的流程一般包括圖像分割、圖像輪廓提取、曲線特征點提取、曲線擬合等主要環(huán)節(jié)。為了更好地實現(xiàn)矢量化,針對各個環(huán)節(jié)存在的突出問題,人們研究和提出了各種不同的解決思路,從而產(chǎn)生了很多不同的矢量化算法,主要包括基于輪廓的矢量化算法、基于邊緣的矢量化算法、基于細(xì)化的算法、基于梯度網(wǎng)格的矢量化算法等。
錢曉峰等[2]提出了一種輪廓點列快速矢量化算法,極大地減少了輪廓需存儲的像素點個數(shù),節(jié)省了內(nèi)存空間,并為進(jìn)一步處理如形狀匹配、編碼等提供了基礎(chǔ)。黃雪優(yōu)[3]采用Canny算法和基于閾值的算法提取邊緣,對邊緣進(jìn)行編輯獲取閉合的輪廓,根據(jù)輪廓點矢量變化率判斷輪廓的特征點并完成輪廓分段,對分段輪廓進(jìn)行擬合初步獲得圖像矢量化邊緣。最后采用矢量圖手動編輯的方式對獲得的矢量圖邊緣做出修正,減少特征點個數(shù),提高矢量圖的精度。石文瑩[4]提出了采用基于細(xì)化和自適應(yīng)網(wǎng)格相結(jié)合的矢量化方法,解決了傳統(tǒng)基于細(xì)化方法容易在交叉處產(chǎn)生畸變的情況。楊云等[5]提出了一種線目標(biāo)半自動矢量化的方法。該方法通過在用戶選擇的線目標(biāo)上加一個滑動窗,采用顏色空間轉(zhuǎn)換、k-均值聚類和區(qū)域生長相結(jié)合的方法對窗口內(nèi)的當(dāng)前線目標(biāo)進(jìn)行自適應(yīng)分割;再通過不斷地沿前進(jìn)方向移動窗口和對窗口內(nèi)的線目標(biāo)進(jìn)行分割、細(xì)化及序貫跟蹤,來完成線目標(biāo)的矢量化。王紅巖等[1]針對煤礦地質(zhì)剖面圖提出了采用不同的算法刪除各種巖性符號,只保留各個地層的邊界線,最后自動矢量化的算法。這些算法在實際問題解決中都得到了應(yīng)用,一些矢量化軟件采用了其中的部分算法,但是對于一些特殊問題,這些算法還存在不足。
文中在深入研究已有矢量化算法的基礎(chǔ)上,針對油藏地質(zhì)剖面圖的特點,提出了一種基于顏色聚類結(jié)合改進(jìn)的Canny邊緣檢測的油藏地質(zhì)剖面圖矢量化方法。
根據(jù)矢量化的一般流程,確定地質(zhì)剖面圖矢量化的流程,如圖1所示。
圖1 地質(zhì)剖面圖矢量化流程
油藏地質(zhì)剖面圖的地層層位之間一般采用不同的圖例填充,其中圖例由顏色和圖案組成。對圖像首先進(jìn)行顏色空間轉(zhuǎn)換,然后可以采用顏色聚類消除層位之間的圖例。
1.2.1 RGB轉(zhuǎn)Lab
常見的顏色模式有RGB模式、HIS模式、HSV模式和Lab模式等。由于采用RGB模式時,各種顏色在數(shù)值上的差距不能很好地反應(yīng)其色彩差異,于是先將圖像由RGB模式轉(zhuǎn)換到Lab模式。Lab模式是國際照明委員會于1976年制定的一種色彩模式,是一種均勻的顏色空間,也是經(jīng)常用來描述人眼可見的所有顏色的最完備的色彩模型。Lab顏色空間由L、a和b三個坐標(biāo)軸組成。L表示亮度,值域是從0到100逐漸由黑變白;a和b的取值范圍都是-128~+127,a值從小到大的變化是指從綠色變?yōu)榧t色,b值從小到大的變化是從藍(lán)色變?yōu)辄S色。一般先將RGB轉(zhuǎn)化到XYZ空間,然后再由XYZ空間轉(zhuǎn)化到Lab色彩空間,見式1~3。
(1)
(2)
(3)
設(shè)在Lab色彩空間模式下的兩種顏色分別為C1和C2,則C1與C2的顏色色差定義為在L、a、b三個通道的歐氏距離,即:
ED(C1,C2)=
(4)
1.2.2 確定初始聚類中心
假設(shè)油藏地質(zhì)剖面圖的圖例庫為Library(其中Library是由圖例填充時的顏色值L、a、b組成,L、a、b分別為Lab顏色模型的三個分量),根據(jù)圖例庫Library分析油藏地質(zhì)剖面圖地層輪廓的填充顏色,從而確定初始聚類中心的顏色值。對于給定的油藏地質(zhì)剖面圖Image(RGB轉(zhuǎn)Lab后的圖像)中的每一個像素點p,計算其與每種圖例顏色Library(C)(C為圖例庫中的一種顏色)的歐氏距離,每個像素點取對應(yīng)最小歐式距離的點的顏色值建立顏色量化圖像Image2,即:Image2=Library(C)。其中任意C1≠C,滿足:ED(Library(C),Image(p)) 1.2.3 顏色聚類 已知k個初始顏色聚類中心,循環(huán)執(zhí)行以下兩步:(1)對每個像素點p,計算其與每個初始聚類中心的歐氏距離,并將該像素劃分到距離最近的聚類中心上;得到該像素的類別Classp∈{1,2,…,k},即將像素點p劃分到第Classp類。(2)重新計算聚類中心。新的聚類中心為圖例庫Library中與類別平均色最接近的顏色。對于每個類別,先計算該類所包含的所有像素的平均色,然后計算該平均色與Library中的每一種顏色的距離,將距離最小值設(shè)為新的聚類中心;直到聚類中心未發(fā)生變化或達(dá)到最大迭代次數(shù)。 1.2.4 平 滑 顏色聚類后圖像中還存在一些孤立的顏色塊。采用圖的最小割算法進(jìn)行平滑處理,通過最小化能量函數(shù)[6]實現(xiàn),重新分配各像素的類別。最終得到顏色聚類分割后的油藏地質(zhì)剖面圖。 邊緣檢測是圖像處理和計算機(jī)視覺中的重要組成部分。邊緣檢測的目的是識別包括圖像特征信息和形狀質(zhì)量的大部分的亮度變化。眾所周知的經(jīng)典邊緣檢測算子包括Log、Robert、Prewitt、Laplace、Sobel和Canny[7]等。不同的算子有不同的優(yōu)勢。由于Canny算子具有更好的信噪比和檢測精度,因此,Canny算子成為評估其他方法的基準(zhǔn),其基本步驟如下: (1)利用一維高斯函數(shù),對圖像進(jìn)行低通平滑濾波; (2)用一階偏導(dǎo)有限差分計算平滑后圖像中各點的梯度值和梯度方向; (3)對梯度幅值進(jìn)行非極大值抑制; (4)用雙閾值算法檢測和連接邊緣。 但Canny算子也存在如下缺點: (1)傳統(tǒng)的Canny邊緣檢測利用高斯函數(shù)進(jìn)行平滑濾波時會過度平滑圖像,使圖像變得模糊,可能扭曲和丟失弱邊緣,特別是邊緣的連接和拐角處;同時,高斯函數(shù)的方差需要人為設(shè)定; (2)在計算像素點的梯度幅值時,對噪聲較敏感,容易檢測出假邊緣或丟失一些真實邊緣的細(xì)節(jié)[8]; (3)在高低閾值的選取時,由于高低閾值不是由圖像邊緣的特征信息決定的,而是需要人為事先設(shè)定好的,因此不具有自適應(yīng)性。 針對以上缺點,Cheng Y[9]用非線性擴(kuò)散濾波器代替高斯函數(shù)平滑圖像,求梯度幅值時考慮了像素對角線方向,最后用Otsu算法自動獲取閾值,在邊緣檢測中取得了更好的精度。Zhao M等[10]針對傳統(tǒng)Canny方法使用高斯函數(shù)去噪時,平滑度取決于高斯核的大小且只對于高斯噪聲效果好的缺點,引進(jìn)了DCT壓縮域的概念,提出通過DCT變換將2D離散像素轉(zhuǎn)換為連續(xù)DCT域,得到DCT系數(shù)和估計噪聲,對DCT系數(shù)進(jìn)行修正,然后通過IDCT獲得平滑圖像的方法,保留了更多有效的邊緣信息。Ma X等[11]為解決傳統(tǒng)Canny邊緣檢測算子過度平滑圖像和適應(yīng)性弱的問題,改進(jìn)了參數(shù)Sigma和高低閾值的獲取方法。Biswas R等[12]針對Canny算法對于模糊圖像高低閾值難于確定的缺點,提出了基于type-2模糊集概念的算法來自動確定閾值。尚長春等[13]采用一種新的四階偏微分方程的降噪算法對圖像去噪,進(jìn)一步提高降噪效果;采用自適應(yīng)閾值的方法對圖像邊緣進(jìn)行檢測,實現(xiàn)了雙閾值的自適應(yīng)提取;基于模糊判決的理論,提出了一種有效的邊緣連接方法。溫陽東等[14]提出了一種基于改進(jìn)非極大值抑制(NMS)過程和雙閾值求取方法的Canny邊緣檢測算子,能夠獲得更好的圖像邊緣。 在研究以上改進(jìn)算法的基礎(chǔ)上,結(jié)合油藏地質(zhì)剖面圖矢量化的特點,即對地層水平層位輪廓矢量化,提出了對Canny算法改進(jìn)的思路,步驟如下: 1.利用一維高斯函數(shù),對圖像進(jìn)行低通平滑濾波(采用Sigma為1;Gaussian mask size 為5)。 2.用Sobel算子對圖像進(jìn)行卷積,求像素點的梯度,采用水平、45°和135°方向上的一階梯度模板,如圖2所示。 圖2 一階梯度模板 三個方向上的一階梯度分量分別為Px(i,j)、P45°(i,j)、P135°(i,j),由圖2中對應(yīng)的一階梯度模板與圖像進(jìn)行卷積得到(其中i,j分別表示圖像像素點的x坐標(biāo)和y坐標(biāo))。每個像素點的梯度幅值和梯度方向計算分別為: (5) (6) 但是求完梯度后依然會存在一些噪聲,采用如下方法消除部分噪聲。 (2)Matrix為梯度圖像中每一點的梯度平均值累加和矩陣,其中每一點存放了該點所經(jīng)過的窗口的梯度平均值的累加和,Matrix1為梯度平均值的累加次數(shù),計算Matrix中每一點的平均值為Matrix2。 (3)判斷如果M中某一點的梯度值小于Matrix2中對應(yīng)點梯度值的b倍(b為倍數(shù)因子),則設(shè)該點的梯度值為0,否則該點的梯度值不變。(采用鄰域窗口20行20列,step為5,倍數(shù)因子為2)。 3.對梯度幅值進(jìn)行非極大值抑制。 4.用雙閾值算法檢測和連接邊緣,其中的高低閾值選取方法采用改進(jìn)的Otsu自動閾值選取方法[15],連接邊緣后,刪除邊緣長度小于等于2的邊緣。 經(jīng)過邊緣輪廓提取后的油藏地質(zhì)剖面圖只包含水平的地層層位輪廓,需要將其矢量化,以便于后期的編輯和應(yīng)用。采用Two-pass算法先對輪廓信息進(jìn)行連通區(qū)域標(biāo)記,然后對標(biāo)記的每個區(qū)域采用Douglas-Peucker矢量數(shù)據(jù)壓縮算法提取特征點(文中采用的DP閾值為1),最后采用三次貝塞爾曲線對提取的特征點進(jìn)行擬合,從而完成油藏地質(zhì)剖面圖的矢量化。 為了驗證提出算法的可行性,實驗在MATLAB中以Lena圖像為例,并分別添加椒鹽噪聲(噪聲密度為0.05)、高斯噪聲(均值為0、方差為0.01)以及泊松噪聲(噪聲密度為0.05)來進(jìn)行驗證,并與傳統(tǒng)經(jīng)典算法Log、Canny邊緣檢測算法的檢測效果作比較。從定量的角度分析,分別計算各邊緣檢測算子所對應(yīng)的峰值信噪比(peak signal-to-noise ratio,PSNR),如表1和表2所示。 表1 不同噪聲環(huán)境下各算法對應(yīng)的PSNR dB 表2 不同噪聲密度下各算法對應(yīng)的PSNR dB 上述實驗結(jié)果表明,針對不同類型的噪聲以及不同程度的高斯噪聲,文中算法具有較好的抗噪性和準(zhǔn)確性。 利用提出的油藏地質(zhì)剖面圖矢量化流程和改進(jìn)算法,在MATLAB中用實際地質(zhì)剖面圖進(jìn)行了矢量化實驗,如圖3和圖4所示。 圖3 實驗結(jié)果1 實驗過程和結(jié)果如下: (1)采用傳統(tǒng)的Canny邊緣檢測算法,對地質(zhì)剖面圖進(jìn)行檢測,結(jié)果見圖3(b)和圖3(e)??梢钥闯觯瑘D3(b)噪聲很多,效果很差。采用改進(jìn)的Canny邊緣檢測算法,對地質(zhì)剖面圖進(jìn)行檢測,結(jié)果見圖3(c)和圖3(f)??梢钥闯鲈肼暽?,較好地保留了水平的層位輪廓,效果好。 (2)對原始圖先進(jìn)行顏色聚類,結(jié)果見圖4(a)和圖4(e);然后再采用改進(jìn)的Canny邊緣檢測算法進(jìn)行邊緣檢測,得到輪廓檢測結(jié)果,見圖4(b)和圖4(f);最后再進(jìn)行矢量化,結(jié)果見圖4(c)和圖4(g)。 (3)對圖4(b)和圖4(f)采用美國斯坦福大學(xué)的Vector magic軟件進(jìn)行矢量化,結(jié)果見圖4(d)和圖4(h)。可以看出,文中矢量化結(jié)果與Vector magic基本相同,但是文中矢量化方法速度更快,得到的曲線特征點更少,更適合矢量化結(jié)果的存儲和編輯。與文獻(xiàn)[1]中針對每一種圖例進(jìn)行刪除的算法相比,文中算法實用性更強(qiáng)。 圖4 實驗結(jié)果2 油藏地質(zhì)剖面圖中的巖性復(fù)雜多樣,針對其特點提出了油藏地質(zhì)剖面圖的矢量化流程,增加了顏色聚類分割,改進(jìn)了Canny邊緣檢測算法。實驗結(jié)果表明,該方法很好地解決了油藏地質(zhì)剖面的矢量化問題,矢量化速度快、質(zhì)量高,為在隨鉆地質(zhì)導(dǎo)向綜合決策中應(yīng)用油藏地質(zhì)剖面奠定了基礎(chǔ)。下一步包含斷層的油藏地質(zhì)剖面的矢量化將是研究的重點。1.3 改進(jìn)的自適應(yīng)Canny邊緣檢測算法
1.4 邊緣輪廓矢量化
2 實驗結(jié)果及分析
2.1 實驗1
2.2 實驗2
3 結(jié)束語