基于NSCT域鄰域收縮的SAR圖像去噪

鐘微宇, 沈汀

ZHONG Weiyu1,2, SHEN Ting 2

1.中國科學院研究生院，北京 100049

2.中國科學院對地觀測與數字地球科學中心，北京 100094

1.Graduate University of Chinese Academy of Sciences, Beijing 100049,China

2.Center for Earth Observation And Digital Earth, Chinese Academy of Sciences, Beijing 100094, China

1 引言

合成孔徑雷達(Synthetic Aperture Radar,SAR)是一種主動式微波遙感器, 可以全天時、全天候成像,在軍事和民用方面發揮著越來越巨大的作用。SAR的成像機理致使SAR圖像存在固有的相干斑噪聲，這將非常不利于 SAR圖像中景物的自動解譯，特別是噪聲較強時，圖像中點目標、邊緣輪廓和紋理細節將被噪聲淹沒。因此SAR圖像中的相干斑噪聲的去除是 SAR圖像處理中急需解決的一個關鍵問題。去噪效果的好壞直接影響到SAR圖像的后續處理，如分割、目標檢測等。

近年來，SAR圖像相干斑噪聲抑制技術得到飛速發展，主要分為兩大類:空域濾波和頻域濾波。空域濾波有Lee[1]、Gamma_MAP[2]等，這些濾波方法在一定程度上可以有效地減弱噪聲的影響，但此類濾波方法往往會對圖像產生過平滑作用，使圖像變模糊，細節信息丟失嚴重。頻域濾波有基于小波變換[3]的去噪方法和基于多尺度幾何變換的去噪方法[4]。小波變換具有時頻局部化特性,它可以從不同的分辨率空間來描述圖像的局部特征，但對于具有線狀奇異的目標函數（如圖像的邊緣），小波系數不再稀疏。因此小波變換在SAR圖像去噪中并不能夠很好地保持圖像中的邊緣輪廓細節信息，其根本原因在于小波分析在二維空間并不是最優的函數表示方法，并不能很好地刻畫圖像中具有線奇異的幾何信息。針對小波變換的不足，2006年，Cunhua等人提出了一種新的多尺度幾何分析工具-非下采樣輪廓變換(NSCT)[5]，用于SAR圖像去噪具有明顯優勢。NSCT通過構造輪廓基來檢測小波基所不能充分刻畫的幾何結構，同時很好地保持了圖像的奇異信息，有利于增強去噪結果中細節結構的保持能力。NSCT繼承了小波變換的時頻特性，同時具有良好的多方向性、各向異性和平移不變性。

2005年，A.Krzyzak[6]等人提出了基于小波域鄰域收縮（NeighShrink）的去噪方法。該去噪方法利用小波系數相關性，使用鄰域收縮法對小波系數進行處理，取得較好去噪效果。后來出現各種改進的NeighShrink去噪方法。2008年，周登文[7]等人針對NeighShrink在各個小波子帶上均使用次優的閾值以及固定的鄰域窗口尺問題，運用 Stein的無偏風險估計改進NeighShrink方法。2010年，宮霄霖[8]提出了基于相關度分析的自適應鄰域小波去噪算法，來解決 NeighShrink不能自動選取合適鄰域窗口大小的問題。2011年，武海洋[9]等人則采用靜態小波變換域中NeighShrink方法,以此來克服正交小波變換域使用全局閾值不能反映小波系數隨尺度變化的不足。上述的這些去噪方法都是在小波域中使用 NeighShrink方法，然而小波變換不具有多方向性、各向異性和平移不變性，致使在小波域使用NeighShrink去噪不能較好地保護紋理細節信息。本文將借鑒鄰域收縮的方法，提出一種基于NSCT域鄰域收縮的 SAR圖像去噪方法。該法利用 NSCT系數的相關性，使用 NeighShrink對不同子帶系數進行收縮，以提高去噪性能和較好地保護邊緣信息。真實SAR圖像的仿真實驗證明了本算法的有效性。

2 基于NSCT域鄰域收縮的SAR圖像去噪

2.1 非下采樣輪廓變換（NSCT）的基本原理

非下采樣輪廓變換（NSCT）是一種利用迭代非下采樣濾波器組來實現一系列多分辨率、多方向和平移不變的頻域子圖像的算法。該算法先將圖像進行多尺度分析，然后再進行多方向分析，兩步分開進行逐層迭代。NSCT的結構為非下采樣金字塔濾波器組(Nonsubsampled Pyramid Filter Bank ，NSPFB)和非下采樣方向濾波器組(Nonsubsampled Directinnal Filter Bank，NSDFB)兩部分，其實現結構原理如圖1所示,其中圖1(a)為NSCT的結構圖(黃色表示通帶，灰色表示阻帶)，圖1(b)為NSCT非下采樣塔形分解和非下采樣方向濾波的頻率分割示意圖。其分解過程可表述為：先利用非下采樣塔式濾波器組（NSPFB）對圖像進行多尺度的分解，得到各種不同頻率的高頻子帶圖像和一個低頻子帶圖像，然后再利用非下采樣方向濾波器組（NSDFB）對得到的各高頻子帶圖像進行多方向分解，從而得到不同尺度、不同方向的子帶圖像（系數）。

圖1 實現NSCT的結構原理圖

2.1.1 非下采樣塔式濾波器組（NSPFB）

NSCT的多尺度特性是通過二通道非下采樣濾波器組(NSPFB)分解獲取的。NSPFB能夠達到類似于 LP的子帶分解結構，第l層分解的低通濾波器的理想帶通支撐為[-p/2l,p/2l]2，相應高通濾波器 H1(z)的理想支撐應為[-π/2l-2,π/2l-1]2[-π/2l,π/2l]2，是低通濾波器的補集。NSCT采用的NSPFB為二通道非下采樣濾波器組，如圖2所示, 分解濾波器{H0(z), H1(z)}和合成濾波器{G0(z), G1(z)}滿足Bezout恒等式：

（1）式保證了 NSPFB滿足完全重構（Perfect Reconstruction,PR）條件。其中，NSPFB的實現去除了 LP的下采樣操作，并對濾波器進行相應的插零上采樣來完成的。當NSPFB分解層數為L時，經分解得到 1個低通子帶圖像和 L個高通子帶圖像，由于沒有下采樣，所有子帶圖像與原始圖像尺寸相同。非下采樣塔式濾波器組的三級分解原理圖如圖3所示，將圖像分為一個低頻子圖（Y0）和三個高頻子圖（Y1、Y2、Y3）。

圖2 NSPFB

圖3 非下采樣塔式濾波器組分解圖

2.1.2 非下采樣方向濾波器組（NSDFB）

NSCT采用的 NSDFB也是一組二通道非下采樣濾波器組，來實現頻域方向分解，如圖4所示。其中分解濾波器{U0(z), U1(z)}和合成濾波器{V0(z), V1(z)}也滿足Bezout恒等式：

（2）式保證NSDFB滿足完全重構條件。NSDFB實現過程是先對濾波器U0(z)和U1(z)進行上采樣，再對上一級二通道方向分解后的子帶圖像進行濾波，從而實現頻域中更為精確的方向分解。某尺度子帶圖像經過l級方向分解后，就得到了2l個與原圖像尺寸大小相同的方向子帶圖像。實現四通道方向NSDFB分解圖如圖5所示。

圖4 NSDFB

圖5 扇形濾波器實現四方向分解示意圖

將 NSPFB與 NSDFB相結合，就可以實現NSCT。原圖像經NSPFB分解后得到的帶通子帶圖像，再輸入到 NSDFB中可獲得圖像的帶通方向信息，從而實現對圖像的多尺度、多方向分解。

設輸入圖像為 f(x, y),大小為M*N，NSCT對圖像的分解過程可由式(3)所示：

其中aJ為低頻子帶，bj,k為j尺度,k方向的高頻子帶。一個NSCT低頻系數可用aJ(m, n)表示，一個NSCT高頻系數可用 bj,k(m, n)表示。其中J為NSPFB的分解層次；lj為第j層的NSDFB的分解層次；j表示 NSPFB分解后的尺度標號；k表示NSDFB中方向標號；m、n表示系數在方向子帶中的空間位置。圖6給出了一幅圖像（zoneplate）的NSCT一層二方向的分解示意圖，金字塔分解為1級，方向濾波數為2。

圖6 zoneplate圖像的NSCT例圖

2.2 鄰域收縮法的原理

針對閾值去噪法對于許多可能包含有用細節信息的變換系數的“過扼殺”（如硬閾值）[10],或是對于噪聲“過保留”（SURE閾值）[11]的問題,文獻[6]提出一種新的去噪方法-鄰域收縮法(NeighShrink)。鄰域收縮法的思想是假定在一個變換域中，在一個較小的系數鄰域內，變換系數具有一定的相關性，即在幅值較大的系數周圍存在較大的系數的概率非常大，所以在對噪聲進行閾值處理時，可以參考鄰域系數的情況，利用尺度內系數相關性。在估計無噪聲變換系數時考慮當前系數鄰域系數的分布情況，而不在是孤立地看待每個系數，如圖7所示。

圖7 鄰域收縮系數圖

NeighShrink去噪算法可由以下幾個步驟來描述：

1）.對含噪圖像進行某種變換，如小波或多尺度幾何變換，將圖像分解為不同頻率的子帶圖像。

2）.對于分解得到 j尺度,k方向的高頻子帶系數 bj,k(m, n)，計算出以(m, n)為中心的窗口內所有系數的平方和，即：

其中，W是以（m,n）為中心的窗口;

3）.計算收縮閾值：

其中 M*N為圖像的大小，σj,k表示圖像j尺度,k方向的高頻子帶的噪聲標準差；

4）.確定收縮因子 βj,k(m, n):

5）.進行系數收縮。收縮后的系數為：

6）.最后對修改后的系數進行相應的逆變換，得到去噪后的圖像。

2.3 結合鄰域收縮的NSCT域SAR圖像去噪算法實現

文獻[5]中 NSCT去噪方法采用硬閾值，對NSCT變換的全部系數采用同一個閾值，通過設定閾值，然后比較系數和閾值，小于閾值的系數將其置為 0。但該法把一部分代表紋理細節信息的系數也置零，導致過度扼殺系數，造成了偽影和模糊。文獻[6]在小波域使用 NeighShrink也有局限性，這是由于小波變換的方向有限，不具有平移不變性。對此，本文將利用NSCT的多方向性、各向異性和平移不變特性，在NSCT域使用NeighShrink處理含噪的變換系數。本文對一幅 river圖像(圖 9(p1))經對數變換后，再進行NSCT分解，分解層數為3，每層方向為 4，4，8，對第二、三層的各一個方向的子帶進行直方圖統計，如圖8(a)和圖 8(b)所示。從直方圖中可看出，圖像經NSCT分解后，其系數往往呈現稀疏性分布：大部分變換系數的幅值較小，它們對應圖像中的光滑部分或白噪聲；只有小部分的幅值較大，這些系數含有一些重要的信息，如圖像的邊緣或奇異位置。在圖像的其他子帶中，也可觀察到類似的分布。

圖8 river SAR圖像NSCT的系數統計直方圖

D.D.Po[4]等人根據鄰域的相關熵理論對輪廓系數的尺度間、尺度內和方向間的相關性進行了詳細分析，并得出結論：輪廓系數在尺度間、尺度內和方向間都存在較強的相關性，其中尺度內變換系數的相關性最強；并將輪廓系數尺度間的相關性與小波系數尺度間的相關性進行了比較，得出前者比后者有更強的相關性。文獻[12]依據輪廓系數尺度間的相關性，提出一種基于雙變量閾值的NSCT圖像去噪方法，實驗結果表明，去噪圖在邊緣特征方面保持了良好的視覺效果。本文則利用 SAR圖像經NSCT分解后，其系數在鄰域內有較強的相關性，提出一種基于NSCT域的鄰域收縮去噪方法。根據鄰域窗口內所有NSCT系數的平方和的大小來決定處于該窗口中心的NSCT系數是收縮還是置零。對少量幅值較大的系數，由尺度內的相關性可得其窗口內平方和較大，由公式6對其進行收縮處理；而對大量幅值較小的系數，當窗口內平方和很小時，對窗口中心系數進行置零處理。相對于軟硬值去噪方法而言，該法由于考慮它們鄰域系數的分布情況，這樣能減少重要的變換系數被誤置為零的情況，可以盡量多地保留圖像細節。

基于 NSCT域鄰域收縮的 SAR圖像去噪算法的步驟由下：

(1).對原始 SAR圖像進行對數變換，把服從Gamma分布的乘性噪聲轉換成加性高斯噪聲。

(2).按照實驗中的 NSCT分解層數和方向分解層數，對經對數變換后的 SAR圖像進行 NSCT分解，得到低頻子帶和高頻子帶圖像（系數）；

(3).由于低頻子帶圖像aJ含有部分噪聲，為更好地平滑噪聲，將用中值濾波去除這部分噪聲，設去噪后的低頻子帶圖像為；

(4).使用魯棒的中值估計子[13]估計不同子帶系數的噪聲標準差σj,k，即

再由式(5)計算 j尺度,k方向的高頻子帶收縮閾值Tj,k。其中Median表示取中值運算，定義為對變換系數按大小排列，再取中位數；

(5).使用式(6)和式(7)對高頻子帶系數進行收縮去噪處理，以估計無噪聲高頻系數；

3 實驗及分析

為驗證本文去噪算法的有效性，將選取邊緣等細節信息比較豐富，具有一定代表性的真實 SAR圖像進行仿真實驗。本文采用二幅SAR圖像來進行驗證，如圖9所示。使用不同去噪算法對這兩幅圖像進行去噪，最后對去噪結果作對比。兩幅圖像的灰度級都是 256，其中第一幅圖像（river）大小為376×376像素，該成像區域主要有河流和陸地這兩部分；第二幅圖像(horsetrack)大小為384×384(美國新墨西哥地區，ku波段1M分辨率)，該區域主要有農作物、草地和光滑地表，這幅圖像可從網站[14]下載。將本去噪算法與小波軟閾值去噪算法[15]、小波鄰域收縮去噪算法[10]、NSCT硬閾值去噪算法[5]用于同一幅SAR圖像去噪，再比較SAR圖像去噪性能量化評估指標。SAR圖像去噪性能量化評估指標[16]采用均值(MEAN)、標準差(STD)、等效視數（ENL）、均值保持指數（PM）、相干斑抑制指數（F）這五個參數，來判斷去噪算法的優劣。

圖9 實驗的二幅原始SAR圖像

在本文去噪算法的實驗參數選擇方面，塔式濾波器組選用‘dmaxflat7'，方向濾波器組選用‘maxflat'；塔式濾波器組采用三層分解[17]，方向濾波器組的分解層數為[2 3 3]。中值濾波器[18]的常用窗口大小有5×5和3×3，經實驗分析，其窗口大小選為 5×5，用其對分解得到的低頻子帶進行濾波。而在鄰域收縮法中的窗口大小選擇方面，文獻[6]實驗結果表明 3×3的窗口大小是最優的，當窗口變大時，去噪能力變差，當窗口大小非常小的，去噪能力也不高。因此本文采用 3×3的窗口，這樣在噪聲抑制和保留紋理細節這兩方面能獲得較好的折衷。圖10和圖11為對應于圖9(p1)和圖9(p2)的不同去噪算法的去噪結果圖。

圖10 river圖像去噪結果圖

圖11 horsetrack圖像去噪結果圖

在主觀評價方面，從圖10（a）和圖11（a）中可明顯看出，采用小波軟閾值去噪算法后，圖像中的相干斑受到一定程度的抑制，但圖像邊緣輪廓被嚴重模糊，線性特征與相干斑被一起平滑，損失較多的紋理信息，同時去噪后的圖像均勻區域出現了一些的吉布斯現象，該法的去噪能力有限。圖 10（b）和圖11（b）顯示小波鄰域收縮去噪法對 SAR圖像的均勻區域進行了比較好平滑，但圖像整體存在模糊失真的現象，圖像細節信息丟失比較嚴重,這是由于小波不是二維圖像的最優表示。圖 10(c)和圖 11(c)表明 NSCT硬閾值去噪算法能較好地去除噪聲，但存在比較明顯的劃痕現象，部分邊緣細節也有一定程度的失真。這是由于該法對不同的變換系數采用同一閾值，過度扼殺了一部分代表紋理細節的系數造成的。從以上兩圖可看出去噪效果還不夠明顯，在紋理細節保持方面還有待改進。圖10(d)和圖 11(d)為基于 NSCT域鄰域收縮的去噪算法的去噪結果圖，平滑區域的噪聲去除與以上幾種去噪算法相比較為徹底，邊緣輪廓清晰銳利，連貫清晰，紋理、邊緣（河流線條（圖10(d)）、區域分界線（圖11(d)））和點目標等細節信息得到較好的保護。本文的去噪算法在圖像質量和視覺效果方面相比于以上幾種方法均有所改善，提高了圖像的平滑性和清晰度。該法良好的去噪效果得益于NSCT具有多分辨率、多方向性的特性和NSCT系數之間的相關性。總之，從圖10和圖片11可知，在NSCT域中使用鄰域收縮法對變換系數進行收縮，能較好地去除SAR圖像噪聲和保護紋理細節，這將非常有利于后續的圖像處理。

在客觀評價方面，本文使用等效視數、相干斑抑制指數、均值、標準差、均值保持指數這5個量化指標來評價去噪結果的優劣。表1和表2給出了原始圖像及4種去噪算法對SAR圖像去噪的量化指標。表1中所得的數據的計算區域（較均勻區域塊）大小為90*80，即圖9(p1)中黃色標記區域；表2中的計算區域大小為80*60，即圖9(p2)中黃色標記區域。

表1 river圖像去噪的量化指標

表2 horsetrack 圖像去噪的量化指標

在表1和表2中，從均值保持指數可以看出，這4種濾波方法PM值都非常接近于1，即上述4種去噪算法的均值和原圖像的均值相差不多，基本上能保持原始圖像的灰度強度。這些方法都有良好的圖像均值保持能力，各算法差別不大。均值保持好說明圖像沒有出現光譜信息的畸變，這將有利于圖像的后續處理和解釋。而在STD、ENL、F這三項指標中，本文的基于NSCT域鄰域收縮的去噪算法在這三項指標中有明顯優勢，其 STD比其它 3種去噪算法的STD小，表明其去噪能力最強；4種去噪算法所得的ENL值均高于原始圖像，而用本文算法所得的ENL都比其它3種去噪算法的ENL大（ENL越大，表明圖像上的相干斑越弱，可解譯性越好）；同時用本文去噪算法得到的 F最大，這表明本文的方法能更好地平滑圖像。綜合上述量化評判指標和去噪后的圖像對比，本文提出的基于NSCT域鄰域收縮的去噪算法在去除SAR圖像噪聲和保護紋理細節方面是有效的。

4 結束語

本文實現了一種基于非下采樣輪廓變換（NSCT）域鄰域收縮的SAR圖像去噪算法，算法首先對受相干斑干擾的SAR圖像進行NSCT分解，然后對分解得到的低頻子帶圖采用中值濾波算法，以去除圖像中的低頻噪聲。接著再利用NSCT分解的高頻系數之間有較強的相關性的特征，使用鄰域收縮法對高頻子帶圖的系數進行收縮，對不同系數使用不同的收縮因子，以提高去噪性能指標和保護紋理細節信息能力。真實SAR圖像的去噪實驗結果表明，本文的去噪算法能使SAR圖像得到較好地平滑，同時兼顧保護紋理細節，使得邊緣清晰光滑。該去噪方法具有較好的抗噪性能和良好的邊緣細節保留能力，用于SAR圖像去噪是有效的，在去噪性能和視覺效果方面相對于上述對比方法均有所提高，為后續的SAR圖像處理提供了有效保證。

[1]Lee J S.Refined filtering of image noise using local statistics[J].Computer,Graphics and Image Processing,1981,15(3):380-389.

[2]Barald A, Parmigiani F.A refined Gamma MAP SAR speckle filter with improved geometrical adaptivity[J].IEEE Transactions on Geoscience and Remote Sensing,1995,33(5):1245-1257.

[3]Xie H,Pierce L E and Ulady F T.SAR speckle reduction using wavelet denoising and markov random field modeling[J].IEEE Transactions on Geosciences and Remote Sensing, 2002, 40(10):2196-2212.

[4]Po D D,Do M N.Directional multiscale modeling of images using the Contourlet transform[J].IEEE Transactions on Image Processing, 2006,15(6):1610-1620.

[5]Cunha A L,Zhou J,Do M N.The Nonsubsampled Contourlet Transform:Theory Design,and Applications[J].Transactions on Image Process ing,2006,15(10):3089-3101.

[6]Chen G Y,Bui T D,Krzyzak A.Image denoising using neighboring wavelet coefficients[J].Integrated Computer-Aided Engineering, 2005,12(1):99-107.

[7]周登文,申曉留.鄰域小波系數自適應的圖像降噪[J].中國圖象圖形學報,2008,13(11):2112-2116.

[8]宮霄霖.基于小波變換的不規則鄰域的數字圖像去噪算法研究[D].天津:天津大學，2010.

[9]武海洋,王慧,樊菊.基于靜態小波變換改進鄰域系數法的遙感圖像濾波[J].海洋測繪,2011,31(2):24-27.

[10]Donoho D L, Johnstone I M ,Kerkyacharian G.Wavelet Shrinkage:Asymptopia[J].Journal of Royal Statistics Society,1995,57(4):301-369.

[11]Donoho D L, Johnstone I M.Adapting to Unknown Smoothness Via Wavele Shrinkage[J].J.Acoust.Soc.Am,1995,90(3):1200-1223.

[12]賈建,焦李成,項海林.基于雙變量閾值的非下采樣Contourlet變換圖像去噪[J].電子與信息學報,2009,31(3):532-535.

[13]Chang S G，Bin Y,Varttereli M.Adaptive wavele Thresholding for image denoising and compression[J].IEEE Transaction on Image Processing，2000,9(9),1532-1546.

[14]Horse track near albuquerque,new mexic [EB/OL].[2012-06-8].http://www.sandia.gov/radar/images/horsetrack.jpg.

[15]Donoho D L.Denoising by soft thresholding[J].IEEE Transactions on Information Theory, 1995,41(3):613-627.

[16]王曉軍,孫洪,管鮑.SAR圖像相干斑抑制濾波性能評價[J].系統工程與電子技術,2004,26(9):1165-1171.

[17]楊曉慧,焦李成.非下采樣Contourlet域GCV準則SAR圖像去噪[J].計算機應用研究,2009,26(9):3542-3544.

[18]Tukey,J.W.Exploratory Data Analysis.Reading, MA.Addison-Wesley,1977.