利用信息熵的巖心圖像自適應(yīng)壓縮感知重構(gòu)①

唐國維,劉彥彤,張巖(東北石油大學(xué) 計算機與信息技術(shù)學(xué)院,大慶 163318)

?

利用信息熵的巖心圖像自適應(yīng)壓縮感知重構(gòu)①

唐國維,劉彥彤,張巖
(東北石油大學(xué) 計算機與信息技術(shù)學(xué)院,大慶 163318)

摘 要:針對BCS-SPL算法對巖心圖像進行壓縮感知重構(gòu)的細節(jié)模糊的問題,提出一種利用信息熵的巖心圖像BCS-SPL壓縮感知重構(gòu)算法.采用小波變換對巖心圖像進行稀疏表示,對各子帶進行多尺度分塊,依據(jù)信息熵的大小自適應(yīng)分配采樣率并確定觀測矩陣,通過維納濾波結(jié)合Landweber迭代操作實現(xiàn)重構(gòu).實驗結(jié)果表明,在相同采樣率下,與原始的BCS-SPL算法相比,該算法的重構(gòu)質(zhì)量提高了2-4 dB.

關(guān)鍵詞:巖心圖像; 壓縮感知; 重構(gòu); 自適應(yīng); 信息熵

巖心是油氣田勘探開發(fā)中重要的基礎(chǔ)地質(zhì)資料,在推斷沉積環(huán)境和生儲蓋組合研究中具有不可替代的作用.將巖心樣本通過掃描方式以數(shù)字圖像形式存儲,已成為數(shù)字化油田建設(shè)的重要組成部分.由于多年的累積和不斷新取心,導(dǎo)致巖心數(shù)據(jù)量及其龐大,因此必須對巖心圖像進行壓縮處理.通過對大量典型巖心圖像分析,發(fā)現(xiàn)巖心圖像普遍具有紋理信息豐富的特點,并且對比度很弱.因此,對于巖心圖像的壓縮與普通自然圖像壓縮有著不同的要求[1].目前,國內(nèi)外各類巖心圖像壓縮重構(gòu)算法均基于Shannon/Nyquist釆樣理論[2,3].由于理論框架的原因,基于Shannon/Nyquist采樣定理對巖心圖像進行壓縮重構(gòu),其采樣數(shù)據(jù)具有非常大的冗余性,需要耗費大量的處理時間和存儲空間,且壓縮效果難以保證.

Donoho、Candès及Tao等人建立的壓縮感知[4,5](Compressed Sensing,CS)理論指出,只要信號是稀疏的或者在某一變換空間是稀疏的或可壓縮的,以遠低于奈奎斯特采樣率的速率隨機采樣,仍能夠精確地重構(gòu)原始信號.可見,CS突破了傳統(tǒng)信源編碼架構(gòu)已經(jīng)接近的理論極限,將其用于圖像壓縮可能獲得意想不到的效果.由于直接使用CS方法重整幅圖像的計算量相當(dāng)巨大,Gan將圖像分塊技術(shù)運用到圖像中來,即分塊壓縮感知(Block Compressed Sensing,BCS)[6]方法.BCS方法把整幅圖像分成等尺寸的塊,獨立地對每個圖像塊進行觀測和重構(gòu),這樣大大降低了存儲和計算成本,但是會在低碼率下帶來塊效應(yīng).為此,Sungkwang等人提出BCS-SPL(Block Compressed Sensing-Smooth Projected Landweber)[7]算法,該算法通過高斯隨機矩陣實現(xiàn)采樣,圖像的重構(gòu)策略使用維納濾波結(jié)合Landweber迭代算法實現(xiàn),雖然改善了塊效應(yīng),但是在一定程度上降低了重構(gòu)質(zhì)量.

根據(jù)油田的實際需要,本文將BCS-SPL算法用于解決巖心圖像壓縮與重構(gòu)問題.通過大量實驗發(fā)現(xiàn),其效果并不理想.分析其原因在于巖心圖像普遍包含豐富的目標(biāo)紋理信息,直接使用通用的BCS-SPL算法必然導(dǎo)致重構(gòu)圖像細節(jié)模糊.所以,本文在BCS-SPL框架的基礎(chǔ)上,提出一種利用信息熵的巖心圖像BCS-SPL壓縮感知重構(gòu)算法.該算法在DWT(Discrete Wavelet Transform)域內(nèi),將變換分解后得到的每一級低頻和高頻子帶進行分塊,然后根據(jù)各個子帶的信息熵得到其自適應(yīng)采樣子率進行自適應(yīng)采樣,再通過維納濾波結(jié)合Landweber迭代實現(xiàn)圖像的重構(gòu),達到進一步提高巖心圖像的重構(gòu)質(zhì)量和改善視覺效果的目的.

1　BCS-SPL算法

根據(jù)壓縮感知理論,假設(shè)x為從M個采樣信號y中獲得的長度為N的信號,且M?N,那么,可以從(1)式中恢復(fù)信號x .

其中,x∈RN,y∈RM,即x是一個N維向量,y是一個M維向量,Φ是一個具有采樣率為S= M N的M′ N維的觀測矩陣.由于x的數(shù)量遠大于觀測值y,理論上通過y恢復(fù)x是不可能的.然而,如果x足夠稀疏,就能夠使精確重構(gòu)成為可能[8].

由于圖像數(shù)據(jù)的多維性,采樣過程的維數(shù)N會隨著圖像x的增大而迅速增加,導(dǎo)致存儲觀測矩陣Φ需要巨大內(nèi)存空間,并且重構(gòu)過程會產(chǎn)生相當(dāng)大的計算量.為此,Gan提出分塊CS方法,文獻[6]給出了一個二維圖像的CS范式.在這個技術(shù)中,圖像采樣是通過應(yīng)用塊到塊基的隨機矩陣實現(xiàn)的,重構(gòu)是Landweber迭代結(jié)合平滑操作實現(xiàn)的.由于分塊CS采樣和平滑Landweber迭代重構(gòu)相結(jié)合,所以稱之為BCS-SPL技術(shù).

假設(shè)一幅大小為N′ N的圖像x被分成大小為B′ B的圖像塊,第i個圖像塊的向量表示記為xi,使用觀測矩陣ΦB進行采樣,得到觀測值.

其中,B的大小根據(jù)圖像重構(gòu)的速率和重構(gòu)的質(zhì)量要求綜合決定: 當(dāng)B較小的時候,內(nèi)存占用少且計算速度快; 當(dāng)B較大的時候,圖像的重構(gòu)效果比較好.i =1Kn ,n =N2B2,ΦB是大小為MB×B2的正交觀測矩陣,M=(M×B2)N2,M為對整幅圖像的觀測采樣數(shù).

在文獻[6]中,維納濾波被納入基本Landweber迭代框架中,目的是為了去除塊效應(yīng).本質(zhì)上,這個操作對于Landweber迭代來說,除了固有稀疏性還能夠起到平滑作用.具體地,維納濾波步驟被插入到公式(3)、(4)中的Landweber迭代,可以看出,第i+1次的迭代圖像x(i+1)近似值是通過x(i)得到的.

BCS-SPL算法將圖像分塊進行觀測采樣,從根本上減少了觀測矩陣的存儲量,使重構(gòu)圖像的效率明顯提高,但對圖像分塊進行觀測采樣會割裂圖像的整體信息,而其投影迭代過程中會產(chǎn)生塊效應(yīng),去除塊效應(yīng)會帶來額外的資源消耗和圖像信息丟失,該算法采用維納濾波來去除塊效應(yīng),導(dǎo)致丟失圖像的邊緣和細節(jié)信息.在BCS-SPL算法的基礎(chǔ)上提出的MS-BCS -SPL算法[9]在DWT域內(nèi),對變換分解后的每一級的子帶進行分塊采樣,再通過平滑迭代重構(gòu)圖像.該算法兼顧了CS的計算開銷和圖像的重構(gòu)質(zhì)量,但是直接將其應(yīng)用到巖心圖像壓縮感知重構(gòu)的應(yīng)用中,重構(gòu)的巖心圖像的細節(jié)信息仍然模糊.

2　巖心圖像壓縮感知采樣與重構(gòu)

2.1巖心圖像頻譜分析與信息熵計算

圖像可以看作是一個離散的二維函數(shù),其頻譜│F(U,V)│可以由該圖像矩陣作二維離散傅里葉變換得到.而二維離散傅里葉變換能夠描述圖像紋理近似周期模式的分布規(guī)律[10],所以本文使用基于傅里葉變換的頻譜分析對巖心圖像的頻譜進行分析.圖1(a)和(b)分別給出了巖心圖像及其傅里葉頻譜.

從圖1(b)可以看出,在巖心圖像的傅里葉頻譜中,相對于遠離坐標(biāo)原點的邊緣位置,其靠近坐標(biāo)原點的中間部分并不是特別亮,也就是說,巖心圖像在低頻部分并沒有顯著的優(yōu)勢,而是在高頻部分同樣占有較重要的比例.這說明巖心圖像在經(jīng)過小波變換后,除了低頻子帶,每一級高頻子帶必然會包含更多的、不同方向的重要信息.

圖1　巖心圖像的傅里葉頻譜

由于巖心圖像的高頻成分占有相對較多的比例,那么圖像分塊后,不同子塊紋理不同,即包含信息量不同.由Mallat塔式小波分解理論可知,圖像經(jīng)分解層數(shù)為L的二維離散小波變換以后,分為3L+ 1個子頻帶,即1個低頻子帶和3個高頻子帶.理想的采樣方法是信息量少的塊少采樣,信息量多的塊多采樣,在總采樣率不變的情況下,將有限的資源有效地分配給紋理相對復(fù)雜的圖像塊.因此,本文改進了所有塊都使用相同的采樣率的BCS采樣方法,根據(jù)子塊間紋理結(jié)構(gòu)不同引起信息量的差異的特點,采用信息熵作為紋理信息的度量.通過計算經(jīng)DWT變換后各級高頻子帶的信息熵,得到自適應(yīng)采樣率.

信息熵反映了巖心圖像中平均信息量的多少,即表示了巖心圖像中灰度分布的聚集特征所包含的信息量,將巖心圖像灰度值進行數(shù)學(xué)統(tǒng)計,便可得到每個灰度值出現(xiàn)的次數(shù)及概率.一般情況下,信息熵值越大表明圖像信息保留的程度越好,其攜帶信息量的能力越強[11].定義巖心圖像信息熵的計算公式為:

其中,pi表示巖心圖像的概率密度函數(shù),可利用直方圖近似計算.利用信息熵的計算公式便可計算出巖心圖像的信息熵.

2.2基于信息熵的自適應(yīng)采樣

在本文的改進算法中,把觀測矩陣Φ分成兩個部分: 一個是DWT多尺度變換矩陣Ω,而另一個是多尺度分塊自適應(yīng)觀測矩陣Φ″,即Φ=Φ″Ω.假設(shè)Ω為L 級DWT分解,那么,Φ″是由3L+1個不同的觀測矩陣組成.這時,被分成大小為Bl×Bl的圖像塊的圖像x在l級的低頻和高頻子帶分別通過自適應(yīng)觀測矩陣Φ″進行采樣.自適應(yīng)采樣的實現(xiàn)步驟如下:

1)計算分解層數(shù)為L的DWT變換的l級的采樣子率Sl:

在l級,根據(jù)其分塊大小Bl使用矩陣Φl進行采樣,會產(chǎn)生采樣子率Sl[12].其中,設(shè)DWT基帶子率S0為全采樣率,即S0= 1.若l級的采樣子率Sl,可以得到公式(6):

這里,Wl為l級的采樣子率Sl的加權(quán)系數(shù).加權(quán)系數(shù)Wl可由下式得到.

那么,整個圖像的采樣率為:

由此可知,當(dāng)已知圖像的目標(biāo)采樣率S和加權(quán)系數(shù)Wl后,由(8)式很容易求出S' ,再通過(6)式得到l級的采樣子率Sl.表1給出了在不同的目標(biāo)采樣率S下,分解層數(shù)L=3的DWT變換實現(xiàn)的l級的采樣子率Sl統(tǒng)計;

表1　L3=級DWT變換實現(xiàn)的采樣子率統(tǒng)計

2)計算l級h ,v ,d子帶的信息熵Hlh,Hlv和Hld;

3)計算l級h ,v ,d子帶的自適應(yīng)采樣子率Sla,見公式(9).

根據(jù)自適應(yīng)采樣的實現(xiàn)步驟可以計算出自適應(yīng)采樣子率.表2給出了巖心圖像1在不同的目標(biāo)采樣率S 下,分解層數(shù)L=3時的DWT變換實現(xiàn)的l級h、v和d子帶的自適應(yīng)采樣子率Sla統(tǒng)計.

表2　L3=級DWT變換實現(xiàn)的高頻子帶自適應(yīng)采樣子率統(tǒng)計

根據(jù)表2的l級的h、v和d子帶的自適應(yīng)采樣子率Sla,可以計算出相應(yīng)方向的高頻子帶的觀測采樣數(shù)M .設(shè)MBl表示塊大小為Bl時觀測采樣數(shù),由計算.表3給出了巖心圖像1在不同的目標(biāo)采樣率S下,L=3級DWT變換實現(xiàn)的每一級的h、v 和d子帶的自適應(yīng)觀測采樣數(shù)MBl統(tǒng)計.

表3　L3=級DWT變換實現(xiàn)的高頻子帶的自適應(yīng)觀測采樣數(shù)統(tǒng)計

根據(jù)表3的每一級的h、v和d子帶的自適應(yīng)觀測采樣數(shù)MBl,可以看出觀測采樣數(shù)MBl會根據(jù)分解級數(shù)不同和高頻子帶方向不同而自適應(yīng)改變,充分體現(xiàn)了其根據(jù)圖像塊所包含信息量不同觀測采樣數(shù)不等的自適應(yīng)性,同樣的工作量卻能夠保留更多的邊緣和細節(jié)信息,從而提高巖心圖像的重構(gòu)質(zhì)量.

2.3多尺度BCS自適應(yīng)重構(gòu)算法

在圖像 DWT 稀疏變換域內(nèi),結(jié)合圖像邊緣的3′ 3維納濾波和稀疏提升閾值處理[13]實現(xiàn)重構(gòu).維納濾波在空間域?qū)崿F(xiàn),而平滑和閾值操作在變換域進行.該算法能夠?qū)崿F(xiàn)圖像的快速重構(gòu),在DWT分解的每一級中的每個塊都使用自適應(yīng)觀測矩陣Φ″和Landweber迭代操作.重構(gòu)算法具體步驟如下:

步驟1: 利用MMSE(Minimum Mean Square Error)估計得到xi的近似解xi,從而得到圖像的初始解xi;

1)用3′ 3鄰域的自適應(yīng)維納濾波器去除圖像分塊所帶來的塊效應(yīng);

2)將濾波后的圖像投影在凸集上,可由下式得到:

當(dāng)Φ是正交矩陣的時候,即ΦΦT=1,上式可簡化為:

3)用小波域雙變量閾值[14]對投影結(jié)果進行濾波;

4)將濾波后的圖像再次投影到凸集上;

5)判斷并終止迭代,直到得到最優(yōu)解.

3　實驗仿真與結(jié)果分析

實驗中用到的測試圖像是512×512的兩幅巖心圖像,對其進行利用信息熵的分塊自適應(yīng)采樣與重構(gòu),并與BCS-SPL和MS-BCS-SPL(Multiscale Block Compressed Sensing with Smoothed Projected Landweber)算法進行比較.本文算法、BCS-SPL和MS-BCS-SPL算法均使用雙樹復(fù)小波變換(Dual-Tree Complex Wavelet Transform,DTCWT)[15]作為稀疏基,采樣時使用97雙正交3級DWT作為多尺度變換矩陣Ω.Ω進行l(wèi)級分解時,使用大小為Bl×Bl的圖像塊采樣.該采樣過程使用隨機DCT(Discrete Cosine Transform)SRM(Site Recovery Manager)觀測矩陣[16]實現(xiàn).所有實驗都在MATLAB R2013b環(huán)境下完成.本文算法中,當(dāng)l=1,2,3時,圖像塊的大小分別為Bl=16,32,64,l級的每個方向的高頻子帶的采樣子率Sl都使用表2的計算結(jié)果,然后根據(jù)表3的l級的每個方向的高頻子帶的自適應(yīng)觀測采樣數(shù)MBl得到自適應(yīng)觀測矩陣Φl″.MS-BCS-SPL算法中,當(dāng)l=1,2,3時,塊的大小分別為Bl=16,32,64,每一級的采樣子率都使用表1的計算結(jié)果.BCS-SPL算法中,B=32.圖2、3給出了三種算法重構(gòu)的巖心圖像的部分實驗結(jié)果.由圖可見,當(dāng)采樣子率S=0.2時,本文提出的利用信息熵的分塊自適應(yīng)采樣和多尺度重構(gòu)的圖像質(zhì)量優(yōu)于BCS-SPL算法約4 dB,也優(yōu)于MS-BCS-SPL算法約2 dB.表3給出三種算法對巖心圖像1和巖心圖像2重構(gòu)結(jié)果的峰值信噪比.

圖2　巖心圖像1的3種重構(gòu)算法效果對比(S= 0.2)

圖3　巖心圖像2的3種重構(gòu)算法效果對比(S= 0.2)

表4　三種算法重構(gòu)結(jié)果的峰值信噪比PSNR(dB)

4　結(jié)論

借鑒DWT域的多尺度分塊壓縮感知技術(shù),本文提出一種利用信息熵的巖心圖像BCS-SPL壓縮感知重構(gòu)算法.該算法針對巖心圖像特性,利用DWT的多分辨率和多尺度特性以及信息熵計算使每級分解層上每個方向的高頻子帶上的采樣具有自適應(yīng)性,所以其觀測結(jié)果能夠充分表示巖心圖像的結(jié)構(gòu)特點.因此,使用本文算法對巖心圖像進行壓縮重構(gòu),在重構(gòu)質(zhì)量和視覺效果方面都有所提升.

參考文獻

1Zhan X,Zhang R,Yin D,Huo C.SAR image compression using multiscale dictionary learning and sparse representation.IEEE Geoscience and Remote Sensing Letters,2013,10(5): 1090–1094.

2Xiao C.Reconstruction of bandlimited signal with lost samples at its Nyquist rate—the solution to a nonuniform sampling problem.IEEE Trans.on Signal Processing,1995,43(4): 1008–1009.

3Yang F,Hu J,Li SQ.A total least squares reconstruction algorithm of UWB signals based on sub-nyquist sampling.Journal of Electronics and Information Technology,2010,32(6): 1418–1422.

4Donoho DL.Compressed sensing.IEEE Trans.on Information Theory,2006,52(4): 1289–1306.

5Candès E,Tao T.Near-optimal signal recovery from random projections: Universal encoding strategies.IEEE Trans.on Information Theory,2006,52(12): 5406–5425.

6Lu G.Block compressed sensing of natural images.2007 15th International Conference on Digital Signal Processing.Cardiff.IEEE.2007.403–406.

7Mun S,Fowler JE.Block compressed sensing of images using directional transforms.2009 16th International Conference on Image Processing.Cairo.IEEE.2009.3021– 3024.

8Candès E,Romberg J,Tao T.Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information.IEEE Trans.on Information Theory,2006,52(2): 489–509.

9Fowler JE,Mun S,Tramel EW.Multiscale block compressed sensing with smoothed projected landweber reconstruction.2011 19th European Signal Processing Conference.Barcelona.IEEE.2011.564–568.

10Chen J,Deng M,Xiao PF,Yang MH,Mei XM,Liu HM.Optimal spatial scale choosing for high resolution imagery based on texture frequency analysis.Journal of Remote Sensing,2011,15(3): 492–511.

11徐潔,張國海.基于熵的圖書館網(wǎng)絡(luò)信息安全的風(fēng)險評估.科技情報開發(fā)與經(jīng)濟,2008,18(14):3–5.

12Candès EJ,Wakin MB.An introduction to compressive sampling.IEEE Signal Processing,2008,25(2): 21–30.

13蔣業(yè)文,于昕梅.基于DWT的多尺度分塊變采樣率壓縮感知圖像重構(gòu)算法.中山大學(xué)學(xué)報:自然科學(xué)版,2013,52(3): 30–33.

14Jia J,Jiao LC,Xiang HL.Using bivariate threshold function for image denoising in NSCT domain.Journal of Electronics and Information Technology,2009,31(3): 532–536.

15Selesnick IW,Baraniuk RG,Kingsbury NC.The dual-tree complex wavelet transform.IEEE Signal Processing,2005,22(6): 123–151.

16Do TT,Tran TD,Lu G.Fast compressive sampling with structurally random matrices.IEEE International Conference on Acoustics,Speech and Signal Processing.Las Vegas,NV.IEEE.2008.3369–3372.

Adaptive Compressed Sensing Reconstruction of Core Images Using Information Entropy

TANG Guo-Wei,LIU Yan-Tong,ZHANG Yan
(School of Computer and Information Technology,Northeast Petroleum University,Daqing 163318,China)

Abstract:Aimed at the details vague problem of Block Compressed Sensing-Smooth Projected Landweber compressed sensing reconstruction of core images,a Block Compressed Sensing-Smooth Projected Landweber compressed sensing reconstruction of core images using information entropy is proposed.The method introduces discrete wavelet transform into the sparse representation and conducts multiscale block for each subband,and then adaptively allocates the sampling rates and determines the measurement matrix.The reconstruction can be achieved by Wiener filter combined with Landweber iterative.The experimental results show that the reconstruction quality is improved by 2-4dB compared with that of Block Compressed Sensing-Smooth Projected Landweber under the same sampling rates.

Key words:core images; compressed sensing; reconstruction; adaptive; information entropy

基金項目:①東北石油大學(xué)研究生創(chuàng)新科研項目(YJSCX2015-034NEPU);黑龍江省教育廳科學(xué)技術(shù)研究項目(12521050)

收稿時間:2015-09-24;收到修改稿時間:2015-11-11