胡 曉,朱 磊,崔 琳,潘 楊,劉珂良
(西安工程大學(xué) 電子信息學(xué)院,陜西 西安 710048)
圖像去噪一直是圖像處理領(lǐng)域的經(jīng)典研究課題[1],而相干斑是合成孔徑雷達(SAR,synthetic aperture radar)對分布式目標作相干處理后在圖像上形成的乘性噪聲[2].相干斑噪聲的存在降低了系統(tǒng)對目標的分辨能力,使得SAR圖像的特征檢測、圖像分割和目標分類等解譯工作復(fù)雜化,因此,抑制相干斑噪聲一直是SAR圖像處理的重要步驟[3].
SAR圖像相干斑抑制主要有兩大目標:一是消除圖像中均勻平坦區(qū)域的乘性噪聲;二是盡量保留原有圖像中的邊緣、紋理等結(jié)構(gòu)信息.近幾十年來,國內(nèi)外學(xué)者為了同時滿足這兩大目標,提出了許多經(jīng)典的SAR圖像降斑技術(shù).早期對相干斑噪聲進行處理的多視處理技術(shù)[4],因其大大降低圖像的分辨率而漸漸被淘汰,人們開始更多地研究濾波后處理技術(shù)消除斑點噪聲,大致分為3類:空域濾波[5-7]、各向異性擴散濾波[8-10]和變換域濾波[11-13].空域濾波是最早采用并廣泛用于相干斑的抑制技術(shù),以Lee濾波[5]、Kuan濾波、Frost濾波[6]和其改進[7]為代表,該類算法采用滑動窗口技術(shù)估計窗內(nèi)局部統(tǒng)計特征自適應(yīng)地進行平滑濾波,由于其僅僅考慮了局部鄰域特征而忽略圖像結(jié)構(gòu)信息,無法同時兼顧抑斑性能和邊緣保護.為了克服此缺陷,Buades等[14]于2005年提出針對自然圖像加性白噪聲的非局部均值(NLM,non-local means)算法.該算法采用圖像的塊相似性替代傳統(tǒng)的單像素點相似性來構(gòu)造權(quán)重,從而可以更好地保護邊緣和紋理等特征;文獻[15-18]將非局部平均技術(shù)運用于SAR圖像的相干斑抑制中,取得了不錯的效果;Chen等[16]提出的NL-CV算法,既保護了圖像邊緣又抑制了相干斑,但高灰度值同質(zhì)區(qū)殘留噪聲嚴重;劉書君等[17]提出基于非局部分類處理的SAR圖像降斑方法,首先將圖像分為同質(zhì)區(qū)和異質(zhì)區(qū),然后分別采用加權(quán)平均和3D變換域硬閾值收縮方法進行去噪,達到了較好的抑斑性能和視覺效果;陳世媛等[18]設(shè)計了基于Gabor濾波器的自適應(yīng)均值算法,不僅有效去除了相干斑噪聲,而且很好保護了紋理和邊緣信息.
針對高灰度值同質(zhì)區(qū)殘留噪聲嚴重這一缺陷,文中給出一種基于非局部平均的SAR圖像相干斑抑制算法.由于乘性噪聲使圖像的高灰度值同質(zhì)區(qū)起伏較大,因此文中對原始SAR圖像進行對數(shù)變換和高斯平滑,加快了高灰度值同質(zhì)區(qū)的衰減速度,然后再計算高斯加權(quán)歐式距離作為相似性測量參量;以可分辨原始SAR圖像平坦區(qū)與邊緣區(qū)的變差系數(shù)(CV,coefficient variation)的倒數(shù)作為自適應(yīng)衰減因子,較好地保護了邊緣;聯(lián)合新的相似性測量參量和自適應(yīng)衰減因子構(gòu)造新的加權(quán)系數(shù)對圖像整體信息進行加權(quán)平均濾波.實驗結(jié)果顯示,該算法抑斑的圖像視覺效果更清晰,高灰度值同質(zhì)區(qū)降噪性能和邊緣保護方面都有較大的提高.
假設(shè)一幅含噪數(shù)字圖像表示為V={V(i)|i∈I},矩陣I表示圖像的坐標范圍,對某一像素,估計值NL[V](i)為整個圖像采用非局部平均濾波后所有像素的加權(quán)平均,為
(1)
(2)
與v(Ni)有相似灰度值的所有鄰域像素在加權(quán)平均時被分配較大的權(quán)重,權(quán)重定義為
(3)
式中:參數(shù)h為衰減因子,代表了濾波器的濾波程度;Z(i)是歸一化常數(shù),有
(4)
由于SAR圖像相干斑是乘性噪聲,經(jīng)典NLM方法不能用于SAR圖像去噪,文獻[16]結(jié)合SAR圖像特點提出NL-CV濾波算法,用變差系數(shù)CV代替經(jīng)典NLM算法中的衰減常數(shù),既保護了圖像邊緣也抑制了相干斑,但高灰度值同質(zhì)區(qū)殘留噪聲明顯.為了解決高灰度值同質(zhì)區(qū)去噪困難這一難題,文中給出一種基于非局部平均的SAR圖像相干斑抑制算法.
首先對SAR圖像進行對數(shù)變換和高斯平滑處理,然后計算處理后圖像的高斯加權(quán)歐氏距離作為相似性測量參量,再以可分辨原始SAR圖像起伏程度的CV的倒數(shù)作為自適應(yīng)衰減因子,結(jié)合新相似性測量參量和自適應(yīng)衰減因子求得新權(quán)重,最后利用新權(quán)重對原始SAR圖像進行非局部加權(quán)濾波.
由于SAR圖像固有的噪聲是乘性噪聲,設(shè)一幅SAR圖像在位置i處的觀測值為Y(i),有用回波信號為X(i),相干斑噪聲為S(i),得
Y(i)=X(i)S(i).
(5)
而經(jīng)典的NLM算法是針對加性高斯白噪聲的,因此需首先對服從Gamma分布的乘性噪聲取對數(shù)處理,將乘性噪聲轉(zhuǎn)化為加性噪聲,噪聲分布近似高斯分布[19],得
logaY(i)=logaX(i)+logaS(i).
(6)
式中:a是底數(shù),此處底數(shù)a取10.然后將得到的圖像通過高斯濾波器進一步消除部分噪聲,此時得到的圖像用Gδ{logaY(i)}表示,那么像素i和像素j的相似性用高斯加權(quán)歐氏距離可表示為
(7)
圖1是本文算法與文獻[16]提出的NL-CV濾波算法對簡單仿真圖的相似性測量比較.
(a) 仿真SAR圖像 (b) 本文算法
(c) NL-CV算法圖 1 仿真圖像的相似性測量對比Fig.1 The similarity measurement contrast of the synthetic images
圖1(b)和圖1(c)分別是圖1(a)的A~D窗用本文算法與NL-CV算法的相似性測量比較.由于經(jīng)典NLM算法采用負指數(shù)加權(quán),因此搜索窗內(nèi)像素與中心像素相似程度越高,則其獲得的權(quán)重就越大.圖1(a)中,A、B兩窗的中心像素正好在邊緣線上,所以只有邊緣附近像素與中心像素最相似,則邊緣附近像素權(quán)重較大從而保護了邊緣.從圖1(b)和(c)看出,對數(shù)變換和高斯平滑后邊緣附近的測量值不但明顯低于兩側(cè)同質(zhì)區(qū)且兩側(cè)同質(zhì)區(qū)相對平滑,可以清晰地區(qū)分邊緣和非邊緣,而之前邊緣兩側(cè)同質(zhì)區(qū)起伏較大,灰度值與邊緣區(qū)較為接近,對邊緣保護不利.C、D窗的中心像素分別在高、低灰度同質(zhì)區(qū),比較圖1(b)C窗、D窗和圖1(c)C窗、D窗看出,NL-CV算法C窗的測量值明顯比D窗的高,因此高灰度值同質(zhì)區(qū)去噪不徹底,而本文算法在兩個搜索窗內(nèi)像素的測量值趨于相近,從而保證了不同區(qū)域相近的去噪效果.
傳統(tǒng)NLM算法對圖像進行降噪時,采用常數(shù)h來控制圖像的濾波程度,當h較小時,圖像濾波程度較小,殘留噪聲嚴重;當h較大時,圖像濾波程度較大,容易造成紋理細節(jié)丟失,存在很大的缺陷.CV在空域濾波中可以較好的區(qū)分邊緣區(qū)和同質(zhì)區(qū),因此文中采用CV替代常數(shù)h從而自適應(yīng)地在不同區(qū)域控制圖像的濾波程度.仿真SAR圖像對應(yīng)的CV圖像如圖2所示.
圖 2 圖1(a)的CV圖像Fig.2 The coefficient variation image of Fig.1(a)
由于CV的邊緣區(qū)呈現(xiàn)高灰度而同質(zhì)區(qū)呈現(xiàn)低灰度這一特性,文中定義新的衰減因子為
(8)
式中:α為衰減常數(shù),可調(diào)節(jié)整個算法的濾波強度;CV(i)為方形鄰域內(nèi)觀測值的標準差與均值的比值.
由新的高斯加權(quán)歐式距離與式(8)構(gòu)建的自適應(yīng)衰減因子,組成本文改進算法的加權(quán)系數(shù)為
(9)
(10)
所以文中改進的NLM加權(quán)濾波可以表示為
(11)
圖3是經(jīng)對數(shù)變換和高斯平滑前后仿真圖像的邊緣保護情況,明顯看出,經(jīng)對數(shù)變換和高斯平滑后,圖3(b)邊緣附近像素權(quán)值較大而邊緣兩側(cè)同質(zhì)區(qū)權(quán)值較小,從而更有利于保護邊緣.
(a) 變換前 (b) 變換后圖 3 圖1(a)A窗與B窗經(jīng)對數(shù)變換和高斯平滑前后的加權(quán)系數(shù)對比Fig.3 The weighting coefficient contrast of figure 1 (a) before and after the logarithmic transformation and Gaussian smoothing
對于任意像素點i,選取以像素i為中心的搜索窗像素集合Ω(i),本文算法的具體步驟如下:
(1) 從Ω(i)內(nèi)選取以中心像素i與其他像素j為中心的小尺度相似窗像素集合Y(i)和Y(j);
(2) 依照式(6)和(7),計算Ω(i)內(nèi)像素i和像素j的相似性測量值DGWED(i,j)′;
(3) 計算Ω(i)內(nèi)每個像素的變差系數(shù)并依照式(8)估計衰減因子h′;
(4) 依照式(9)計算Ω(i)內(nèi)每個像素的加權(quán)系數(shù);
(5) 計算Ω(i)內(nèi)各像素灰度值和對應(yīng)加權(quán)系數(shù)的加權(quán)平均值,并作為像素i的濾波估計值.
最后,對整幅圖像的所有像素按上述5個步驟進行濾波,得到濾波后的SAR圖像.
由于非局部平均比局部平均算法計算量大,文中針對某一個像素點進行時間復(fù)雜度分析.設(shè)某一像素點k,其搜索窗尺度為m1×m2,相似窗尺度為n1×n2,則在像素點k的搜索窗內(nèi),在相似性測量部分,由式(6)大約需要n1×n2次乘法,由式(7)大約需要2×n1×n2次乘法和2×n1×n2次加法;計算自適應(yīng)衰減因子部分,由式(8)大約需要n1×n2次乘法;在加權(quán)平均部分,由式(9)和(10)大概需要2×n1×n2次乘法,加權(quán)濾波時,由式(11)需要2×n1×n2次乘法和n1×n2次加法.總體來說,本文算法大約需要8×n1×n2次乘法和3×n1×n2次加法,與NL-CV算法的時間復(fù)雜度基本相似.
為了客觀評價本文算法的有效性和穩(wěn)定性,對圖4所示的仿真圖像和真實圖像分別用本文算法、經(jīng)典Frost濾波[6]和NL-CV濾波[16]進行降噪處理,并從視覺效果和參數(shù)性能2個方面來判別抑斑性能的優(yōu)劣.視覺效果包括去噪圖像和其Canny邊緣檢測圖像,實驗結(jié)果如圖5,6所示.抑斑參數(shù)性能采用了文獻[9]中的方法:等效視數(shù)(ENL,equivalent number of looks)和邊緣保持指數(shù)(EPI,edge preservation index).ENL表示相干斑的平滑能力,其值越大表明相干斑抑制程度越高,本文2次實驗選取的同質(zhì)區(qū)為圖4中的黑色方框;EPI表明圖像邊緣的保護能力,其值越大表明抑斑后圖像的邊緣保護程度越高,理想值為1.不同圖像的EPI和不同區(qū)域的ENL如表1所示,并對最好的性能指標進行了字體加粗.
(a) 3視仿真SAR圖像 (b) 5視真實SAR圖像圖 4 實驗所用SAR圖像Fig.4 Two SAR images for experiment
本文用各算法對2幅圖像進行抑斑時,分別采取了各自最好的抑斑性能參數(shù).各算法對圖4(a)采取的參數(shù)如下:FROST算法的加權(quán)窗尺度為15×15,局部統(tǒng)計量估計窗尺度為7×7,衰減常數(shù)為0.5;NL-CV算法搜索窗尺度為21×21,相似窗尺度為7×7,衰減常數(shù)為300,高斯核標準差為8;本文算法搜索窗尺度為21×21,相似窗尺度為7×7,衰減常數(shù)為0.25;Canny算子的邊緣檢測閾值為0.07.各算法對圖4(b)抑斑時采用的與圖4(a)不同的參數(shù)如下:FROST算法衰減常數(shù)為2;NL-CV算法衰減常數(shù)為100;本文算法衰減常數(shù)為0.06,Canny算子的邊緣檢測閾值為0.06.
FROST NL-CV 本文算法(a) 去噪圖像
FROST NL-CV 本文算法(b) 邊緣檢測圖像圖 5 仿真SAR圖像實驗結(jié)果Fig.5 Experimental results for synthetic SAR image
FROST NL-CV 本文算法(a) 去噪圖像
FROST NL-CV 本文算法(b) 邊緣檢測圖像圖 6 真實SAR圖像實驗結(jié)果Fig.6 Experimental results for real SAR image
對比圖5(a)與圖6(a)可以發(fā)現(xiàn),F(xiàn)ROST濾波的圖像存在模糊,NL-CV算法的抑斑圖像高低灰度同質(zhì)區(qū)殘留噪聲相差較大,而本文算法的抑斑圖像比較清晰,高低灰度值同質(zhì)區(qū)的相干斑噪聲得到了較好地抑制.對比圖5(b)與圖6(b)的邊緣檢測圖像可知,F(xiàn)ROST和NL-CV算法的邊緣檢測圖像在高灰度同質(zhì)區(qū)存在較多的虛假邊緣,而本文算法的邊緣檢測圖像邊緣相對自然平滑,高灰度同質(zhì)區(qū)的虛假邊緣少.
表 1 各算法對2幅SAR圖像抑斑參數(shù)性能對比Table 1 Despeckling parameters contrast of each algorithm for two SAR images
由表1所示的各算法對2幅圖像的抑斑性能參數(shù)比較可以看出,本文改進算法在高灰度值同質(zhì)區(qū)ENL和EPI指數(shù)兩個方面均優(yōu)于其他2種算法,而FROST和NL-CV算法在低灰度同質(zhì)區(qū)的ENL指標更好.FROST濾波因只有局部像素參與運算而抑斑性能較差,NL-CV算法雖采用了非局部思想,但由于乘性噪聲的特殊性使得高灰度值同質(zhì)區(qū)殘留噪聲嚴重,本文算法不僅遵循了非局部思想,而且在相似性測量時通過對數(shù)處理和高斯平滑,有效抑制了高灰度值同質(zhì)區(qū)的噪聲.由于本文算法在相似性測量時對圖像作對數(shù)變換處理,在一定程度上影響了低灰度值區(qū)域的抑斑性能,因為SAR圖像相干斑是乘性噪聲,高灰度值同質(zhì)區(qū)的起伏比低灰度值同質(zhì)區(qū)大,而其他2種算法在低灰度值區(qū)域抑制效果很好,因此加大高灰度值區(qū)域的噪聲抑制是非常有必要的.
針對SAR圖像的乘性相干斑抑制問題,給出了一種基于非局部平均的SAR圖像濾波算法.該算法首先對原始SAR圖像作對數(shù)變換和高斯平滑處理,以加快高灰度區(qū)域的衰減速度從而提高了高灰度同質(zhì)區(qū)的抑斑性能,然后計算處理后圖像每個像素的相似性測量參量.與此同時,采用能較好檢測圖像起伏的CV作為自適應(yīng)衰減因子,從而不影響抑斑性能的同時保護了邊緣.以上兩個步驟使得本文算法的抑斑圖像清晰,提高了高灰度值同質(zhì)區(qū)的相干斑抑制能力和邊緣保護性能.
參考文獻(References):
[1] 楊國梁,雷松澤.基于貝葉斯估計自適應(yīng)軟硬折衷閾值 Curvelet 圖像去噪技術(shù)[J].西安工程大學(xué)學(xué)報,2011,25(6):857-861.
YANG G L,LEI S Z.The image denoising method of soft and hard adaptive thresholding based on Curvelet transform and Bayesian estimation[J].Journal of Xi′an Polytechnic University,2011,25(6):857-861.
[2] WANG Y,AINSWORTH T L,LEE J.On characterizing high-resolution SAR imagery using kernel-based mixture speckle models[J].IEEE Geoscience and Remote Sensing Letters,2015,12(5):968-972.
[3] GAO F,XUE X,SUN J,et al.A SAR image despeckling method based on two-dimensional S transform shrinkage[J].IEEE Transactions on Geoscience & Remote Sensing,2016,54(5):3025-3034.
[4] OLIVER C,QUEGAN S.Understanding synthetic aperture radar images[M].Boston:Artech House,2009.
[5] LEE J S.Digital image enhancement and noise filtering by use of local statistics[J].IEEE Transactions on Pattern Analysis & Machine Intelligence,2009,PAMI-2(2):165-168.
[6] FROST V S,STILES J A,SHANMUGAN K S,et al.A model for radar images and its application to adaptive digital filtering of multiplicative noise[J].IEEE Transactions on Pattern Analysis & Machine Intelligence,1982,4(2):157-165.
[7] 朱磊,水鵬朗,程冬.基于混合迭代濾波的SAR圖像相干斑抑制[J].電子與信息學(xué)報,2012,34(5):1038-1044.
ZHU L,SHUI P L,CHENG D.SAR image despeckling based on mixed iteration filtering[J].Journal of Electronics & Information Technology,2012(5):1038-1044.
[8] YU Y,ACTION S T.Speckle reducing anisotropic diffusion[J].IEEE Transactions on Image Processing A Publication of the IEEE Signal Processing Society,2002,11(11):1260-1270.
[9] ZHU L,ZHAO X T,GU M H.SAR image despeckling using improved detail-preserving anisotropic diffusion[J].Electronics Letters,2014,50(15):1092-1093.
[10] 朱磊,韓天琪,水鵬朗,等.一種抑制合成孔徑雷達圖像相干斑的各向異性擴散濾波方法[J].物理學(xué)報,2014,63(17):445-455.
ZHU L,HAN T Q,SHUI P L,et al.An anisotropic diffusion filtering method for speckle reduction of synthetic aperture radar images[J].Acta Physica Sinica,2014,63(17):445-455.
[11] LI H C,HONG W,WU Y R,et al.Bayesian wavelet shrinkage with heterogeneity-adaptive threshold for SAR image despeckling based on generalized gamma distribution[J].IEEE Transactions on Geoscience and Remote Sensing,2013,51(4):2388-2402.
[12] ARGENTI F,BIANCHI T,LAPINI A,et al.Fast MAP despeckling based on Laplacian—Gaussian modeling of wavelet coefficients[J].IEEE Geoscience & Remote Sensing Letters,2011,9(1):13-17.
[13] PURANIKMATH S S,VANI K.Enhancement of SAR images using curvelet with controlledshrinking technique[J].Remote Sensing Letters,2016,7(1):21-30.
[14] BUADES A,COLL B,MOREL J M.A non-local algorithm for image denoising[J].IEEE Computer Society Conference on Computer Vision & Pattern Recognition.San Diego:IEEE,2005:60-65.
[15] PARRILLI S,PODERICO M,ANGELINO C V,et al.A nonlocal SAR image denoising algorithm based on LLMMSE wavelet shrinkage[J].IEEE Transactions on Geoscience & Remote Sensing,2012,50(2):606-616.
[16] CHEN S B,HOU J H,ZHANG H,et al.De-speckling method based on non-local means and coefficient variation of SAR image[J].Electronics Letters,2014,50(18):1314-1316.
[17] 劉書君,吳國慶,張新征,等.基于非局部分類處理的SAR圖像降斑[J].系統(tǒng)工程與電子技術(shù),2016,38(3):551-556.
LIU S J,WU G Q,ZHANG X Z,et al.SAR image despeckling via the classification-based non-local clustering [J].Systems Engineering and Electronics,2016,38(3):551-556.
[18] 陳世媛,李小將.基于自適應(yīng)非局部均值的SAR圖像相干斑抑制[J].系統(tǒng)工程與電子技術(shù),2017,39(12):2683-2690.
CHEN S Y,LI X J.SAR image despeckling based on adaptive non-local means[J].Systems Engineering and Electronics,2017,39(12):2683-2690.
[19] APRIL G,ARSENAULT H H.Properties of speckle integrated with a finite aperture and logarithmically transformed[J].Journal of the Optical Society of America,1976,66(11):1160-1163.