基于密度聚類與灰度變換的NSST 域聲吶圖像去噪

劉光宇，曾志勇?，曹禹，趙恩銘，邢傳璽

（1.大理大學(xué)工程學(xué)院，云南大理 671003；2.哈爾濱工程大學(xué)物理與光電工程學(xué)院，黑龍江哈爾濱 150001；3.云南民族大學(xué)電氣信息工程學(xué)院，云南昆明 650031）

隨著經(jīng)濟(jì)社會的發(fā)展和相關(guān)技術(shù)的成熟，人們探索海洋的欲望越來越強(qiáng)烈.聲吶探測技術(shù)作為海洋探索的基本手段，在海底地形測繪、水下目標(biāo)物體探測以及生物種群監(jiān)測方面都發(fā)揮著重要作用［1］.由于受海風(fēng)、洋流、水溫、雜質(zhì)以及成像設(shè)備等影響，使用聲吶探測技術(shù)得到的圖像中通常包含各種類型的噪聲，主要表現(xiàn)為顆粒狀的斑點，尤其在淺海區(qū)域，該現(xiàn)象更為嚴(yán)重［2］.斑點噪聲的存在導(dǎo)致圖像質(zhì)量嚴(yán)重下降，對后續(xù)的圖像處理產(chǎn)生不利影響［3］.因此，如何有效去除圖像中的噪聲成了人們研究的熱點.與光學(xué)圖像一樣，聲吶圖像中常見的去噪方法主要包括空間域去噪法（如均值濾波、中值濾波和維納濾波等）和變換域去噪法（如小波變換、脊波變換和輪廓波變換等）［4］.其中變換域去噪法因其復(fù)雜程度低，且具有良好的邊緣保持能力而更受相關(guān)研究者們的關(guān)注［5］.

近年來，剪切波作為變換域中一種新的多尺度幾何分析工具在各種圖像去噪方法中脫穎而出，以其結(jié)構(gòu)簡單，對多維數(shù)據(jù)的各向異性特征近乎最優(yōu)的稀疏表示而深受人們的青睞［6］.由于剪切波變換后圖像存在一定程度的偽吉布斯現(xiàn)象［7］，因此，Easley 等人在其基礎(chǔ)上進(jìn)行改進(jìn)，提出了非下采樣剪切波變換［8］.2020年，Routray等把非下采樣剪切波和雙邊濾波進(jìn)行結(jié)合，構(gòu)建了一種新的圖像去噪模型［9］；同年，Morteza 等提出了基于TLS（t-locations scale）和非下采樣剪切波的噪聲濾波方法［10］.這些方法在保持剪切波原有優(yōu)點的同時，提高了圖像的峰值信噪比和結(jié)構(gòu)相似度，降低了均方誤差，克服了圖像中出現(xiàn)的偽吉布斯現(xiàn)象，使圖像的質(zhì)量得到進(jìn)一步提升.

由于聲吶圖像的噪聲特性與光學(xué)圖像、紅外圖像等存在差異，直接運用光學(xué)圖像和紅外圖像等的去噪方法將難以使聲吶圖像達(dá)到理想的去噪效果.斑點噪聲在聲吶圖像中的分布呈現(xiàn)密度不一的特點，主要集中在混響區(qū)域和目標(biāo)區(qū)域，在陰影區(qū)域密度較?。?1］，與常見的噪聲類型相比，其復(fù)雜性更強(qiáng).聚類算法在圖像處理中已有成功應(yīng)用.其中，DBSCAN（Density-based Spatial Clustering of Applications with Noise）算法可根據(jù)數(shù)據(jù)的分布特性自主計算鄰域值，利用數(shù)據(jù)內(nèi)部的聚集性和連通性可實現(xiàn)更為準(zhǔn)確和高效的數(shù)據(jù)分類［12］.

本文把非下采樣剪切波變換與DBSCAN 密度聚類相結(jié)合，利用灰度變換對聲吶圖像進(jìn)行處理，構(gòu)建一種新的圖像去噪模型.該模型在去除斑點噪聲和緩解偽吉布斯現(xiàn)象的同時，進(jìn)一步保持了圖像的邊緣細(xì)節(jié)特征，對提高圖像峰值信噪比、結(jié)構(gòu)相似度和降低均方誤差等擁有明顯的效果.

1 基本理論

1.1 噪聲模型

根據(jù)噪聲與信號的關(guān)系，圖像的噪聲類型可分為加性噪聲和乘性噪聲［13］.聲吶圖像中常見的噪聲類型為乘性斑點噪聲，其模型可描述為：

式中：f（x，y）為含有噪聲的圖像；n（x，y）為原始圖像；q（x，y）為噪聲信號.乘性斑點噪聲通常服從瑞利分布，與光學(xué)圖像中常見的加性高斯噪聲相比，具有起伏劇烈、均勻度較低以及存在較大不確定性等特點，處理起來更為困難［14］.加性高斯噪聲和乘性斑點噪聲在圖像中的分布如圖1所示.

圖1 噪聲分布模型圖Fig.1 Model diagram of noise distribution

1.2 剪切波變換

剪切波變換（Shearlet Transform，ST）是在合成小波變換基礎(chǔ)上發(fā)展的數(shù)學(xué)算法，可通過對基函數(shù)實行伸縮、平移和剪切等操作生成所需的剪切波函數(shù)，是一種接近最優(yōu)的多維函數(shù)稀疏表示方法［15］.

當(dāng)維數(shù)n=2時，剪切波的仿射系統(tǒng)為：

式中：ψ∈L2（R2）；A、B均為二維可逆矩陣，且|detB|=1；j和l分別為尺度參數(shù)和幾何參數(shù)；k為平移參數(shù).如果AAB（ψ）滿足緊支撐框架，即

則AAB（ψ）系統(tǒng)中的元素ψ就稱為合成小波，其中膨脹矩陣Aj代表尺度變換，矩陣Bl代表幾何變換.通過上述仿射變換，系統(tǒng)可以構(gòu)造各個尺度和各個方向上的緊支撐框架［16］.令式（2）中的A為各向異性的膨脹矩陣，B為剪切矩陣，通常情況下：

當(dāng)a=4，s=1時，

此時的合成小波被稱為剪切波［17］.剪切波是一個帶有尺度、方向和位置的函數(shù)集合.在頻域中，不同尺度的剪切波在其頻域支撐下，表現(xiàn)為近似22j×2j大小的梯形形狀.剪切波基函數(shù)的支撐區(qū)域可以隨著尺度的長寬比發(fā)生變化，能夠最優(yōu)地表示圖像的內(nèi)在幾何結(jié)構(gòu)，其頻域剖分和頻域支撐示意圖如圖2 所示.圖2（a）中的ξ1和ξ2為剪切波頻域剖分的尺度參數(shù).

圖2 剪切波頻域剖分及其頻域支撐示意圖Fig.2 Shearlet frequency domain subdivision diagram and schematic diagram of frequency domain support

1.3 非下采樣剪切波變換

傳統(tǒng)剪切波變換在多尺度分解過程中存在下采樣操作，容易導(dǎo)致細(xì)節(jié)信息丟失，重構(gòu)后的圖像極易產(chǎn)生偽吉布斯現(xiàn)象.非下采樣剪切波變換（Nonsubsampled Shearlet Transform，NSST）的出現(xiàn)解決了上述問題，在保持傳統(tǒng)剪切波良好的局部性、方向敏感性、最佳稀疏性等特性的同時，有效緩解了頻譜混疊現(xiàn)象［18］.

非下采樣剪切波變換的離散化過程分為多尺度分解和方向局部化兩部分［19］：

1）多尺度分解：利用非下采樣金字塔濾波器（Non-subsampled Pyramid Filter，NSPF）將圖像分解為k+1 個與原圖像大小相同的子帶圖像，其中包括1個低頻子帶圖像和k個高頻子帶圖像.

2）方向局部化：標(biāo)準(zhǔn)的剪切波變換中使用的剪切濾波器（Shear Filter，SF）是在偽極化坐標(biāo)系中通過窗函數(shù)的平移操作實現(xiàn)的，此過程中包含下采樣操作，不具備平移不變性.而非下采樣剪切波變換將標(biāo)準(zhǔn)的剪切濾波器從偽極化網(wǎng)絡(luò)系統(tǒng)映射回到笛卡爾坐標(biāo)系統(tǒng)，摒棄下采樣操作，實現(xiàn)了平移不變性.

非下采樣剪切波的具體實現(xiàn)過程：首先，采用非下采樣金字塔濾波器將圖像分解為低頻子帶圖像和高頻子帶圖像.其次，對高頻子帶圖像構(gòu)建Meyer窗，進(jìn)行方向局部化，獲得不同的方向子帶.最后，對每一個方向子帶進(jìn)行傅立葉逆變換，以獲得非下采樣剪切波系數(shù).NSST 二級分解示意圖如圖3所示.

圖3 NSST二級分解示意圖Fig.3 Schematic diagram of NSST level two decomposition

1.4 密度聚類

DBSCAN 密度聚類作為聚類中常見的算法，因其具有識別任意形狀的簇、不需要事先設(shè)置簇的數(shù)量以及能夠識別噪聲點等優(yōu)點，受到人們的青睞［20］.它的基本思想是定義一個簇為所有緊密連接的數(shù)據(jù)點的最大集合，在以某個數(shù)據(jù)點為圓心，以ε為半徑的區(qū)域內(nèi)，DBSCAN 通過區(qū)域查詢，判斷該區(qū)域包含的數(shù)據(jù)點個數(shù)是否大于設(shè)定的密度閾值δ，如果大于該數(shù)值，則創(chuàng)建一個新的簇，然后把相關(guān)聯(lián)的簇進(jìn)行合并.如果某個對象不在任何簇中，則該對象被視為噪聲點，從而達(dá)到分類的目的［21］.基于參數(shù)ε和δ，給出數(shù)據(jù)集D=｛x1，x2，…，xm｝，其中xi∈D，xj∈D.DBSCAN密度聚類算法有以下基本概念［22］：

1）核心對象：若以xj為圓心，ε為半徑的區(qū)域內(nèi)至少包含δ個樣本，則稱xj為核心對象；

2）密度直達(dá)：若xj位于以xi為圓心，ε為半徑的區(qū)域內(nèi)，且xi為核心對象，則稱xj由xi密度直達(dá)；

3）密度可達(dá)：若存在樣本序列K1，K2，…，Kn，其中K1=xi，Kn=xj且Ki+1由Ki密度直達(dá)，則稱xj由xi密度可達(dá)；

4）密度相連：若存在xh均可使xi和xj由xh密度可達(dá)，則稱xi和xj密度相連.

當(dāng)δ=3 時，DBSCAN 密度聚類示意圖如圖4 所示.其中，虛線圈內(nèi)的區(qū)域是以ε為半徑的區(qū)域，x1為核心對象，x2由x1密度直達(dá)，x3由x1密度可達(dá)，x3與x4密度相連.

圖4 DBSCAN密度聚類示意圖Fig.4 Schematic diagram of DBSCAN density clustering

1.5 灰度變換

灰度變換是指將原圖像中的灰度值映射到另一灰度值，其中比例線性變換、分段線性變換是常見的映射函數(shù)［23］.由于分段線性變換在一定的灰度值區(qū)間內(nèi)具有不同的變換函數(shù)，比單一的線性變換函數(shù)更為靈活和更具針對性，可根據(jù)需求拉伸圖像每一部分的灰度值，達(dá)到增強(qiáng)圖像質(zhì)量，提高對比度的效果.分段線性變換函數(shù)示意圖如圖5所示.

圖5 分段線性變換函數(shù)示意圖Fig.5 Schematic diagram of piecewise linear transformation function

2 NSST域聲吶圖像去噪

本文結(jié)合非下采樣剪切波變換、DBSCAN 密度聚類和灰度變換對聲吶圖像進(jìn)行去噪，具體去噪步驟如下：

1）利用非下采樣金字塔濾波器對含噪圖像進(jìn)行4 層分解，分解為1 層低頻子帶圖像和4 層高頻子帶圖像，同時通過剪切濾波器把4 層高頻子帶圖像分別分為23、23、24和24個方向.含噪圖像及其高低頻子帶圖像如圖6所示.

圖6 含噪圖像及其高低頻子帶圖像Fig.6 Noisy image and its high and low frequency subband images

從圖6 中可以看出，含噪圖像經(jīng)過分解后，原始信息被分解到不同層次.其中，低頻子帶圖像包含了原始信息的主要能量，高頻子帶圖像包含了大量細(xì)節(jié)信息和噪聲信號.隨著分解層數(shù)的增加，細(xì)節(jié)信息和噪聲信號越來越集中.因此，將更高層次的子帶圖像分成多個方向，以精準(zhǔn)捕獲每個方向上的細(xì)節(jié)特征.

2）根據(jù)斑點噪聲在聲吶圖像中的分布特性，利用DBSCAN 密度聚類對每一層的高頻子帶圖像進(jìn)行聚類處理，δ的取值通常為N+1，其中N為數(shù)據(jù)的維度.通過調(diào)節(jié)ε值，將信號劃分為細(xì)節(jié)信號和噪聲信號，對噪聲信號進(jìn)行零值處理，保留細(xì)節(jié)信號.

3）低頻子帶圖像中通常不包含噪聲干擾，因此，只需對其對比度進(jìn)行調(diào)整.利用灰度變換中的分段函數(shù)對低頻子帶圖像進(jìn)行灰度拉伸，進(jìn)一步提升圖像的整體視覺效果.

4）利用非下采樣剪切波逆變換對處理后的高頻子帶圖像和低頻子帶圖像進(jìn)行重構(gòu)，得到恢復(fù)后的圖像.

整體去噪流程圖如圖7所示.

圖7 整體去噪流程圖Fig.7 Overall denoising flow chart

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

實驗1選取大小為256×256的沉船圖像［（圖8（a）］作為實驗對象，對其添加方差為0.5 的斑點噪聲，運用本文所提方法對其進(jìn)行處理.得到重構(gòu)后的圖像與變換域中的小波變換（Wavelet Transform，WT）、脊波變換（Ridgelet Transform，RT）、輪廓波變換（Contourlet Transform，CT）、剪切波變換（Shearlet Transform，ST）和非下采樣剪切波變換（Non-subsampled Shearlet Transform，NSST）去噪法的去噪結(jié)果進(jìn)行對比.對比結(jié)果如圖8所示.

圖8 沉船圖像Fig.8 Images of shipwreck

從圖8（c）（d）（e）可以看出，通過小波變換、脊波變換和輪廓波變換去噪后，圖像中的斑點噪聲得到一定程度的去除，但同時圖像整體出現(xiàn)了模糊，細(xì)節(jié)信息丟失嚴(yán)重.在圖8（f）中，通過剪切波變換去噪后，圖像中的噪聲明顯減少，但其邊緣保持能力不足，影響了圖像的質(zhì)量.在圖8（g）中，通過非下采樣剪切波變換處理后，圖像中的偽吉布斯現(xiàn)象得到明顯改善，邊緣保持能力相比剪切波去噪法有了較大的提升.在圖8（h）中，采用DBSCAN 密度聚類、灰度變換與非下采樣剪切波變換相結(jié)合的方法進(jìn)行處理后，圖像中的噪聲得到有效濾除，邊緣輪廓和細(xì)節(jié)特征得到準(zhǔn)確捕捉，失真程度明顯小于小波變換（WT）、脊波變換（RT）、輪廓波變換（CT）、剪切波變換（ST）和非下采樣剪切波變換（NSST）去噪法，從直觀視覺的角度證明了該方法的優(yōu)越性.

為進(jìn)一步客觀具體地評估所提出的方法，對圖8（b）中的含噪圖像（Noisy Image，NI）重復(fù)實驗20次，計算去噪后圖像的質(zhì)量評價指標(biāo)［24］，即均方誤差（Mean Squared Error，MSE）、信噪比（Signal to Noise Ratio，SNR）、峰值信噪比（Peak Signal to Noise Ratio，PSNR）和結(jié)構(gòu)相識度（Structural Similarity，SSIM）的平均值，并通過數(shù)據(jù)對圖像質(zhì)量進(jìn)行判斷.MSE 值、SNR值、PSNR值和SSIM值的計算公式分別如下：

式中：QMSE表示MSE 的值；m×n表示圖像的大?。籌（i，j）為原始圖像；R（i，j）為重構(gòu)后的圖像.MSE值越小表明圖像質(zhì)量越好.

式中：QSNR表示SNR的值.

式中：QPSNR表示PSNR 的值；Lmax表示原始圖像I（i，j）的最大灰度值.與SNR 值一樣，PSNR 值越大表明圖像失真程度越小，圖像質(zhì)量越好.

式中：QSSIM表示SSIM 的值；ux、uy分別為圖像x和y的均值；σx、σy分別為圖像x和y的標(biāo)準(zhǔn)差.SSIM 是衡量兩幅圖像相似度的指標(biāo)，其值越大表明兩幅圖像越相似.

經(jīng)過多次實驗后，沉船圖像的均方誤差（MSE）、信噪比（SNR）、峰值信噪比（PSNR）和結(jié)構(gòu)相識度（SSIM）的平均值如表1所示.

表1 沉船圖像質(zhì)量評價指標(biāo)Tab.1 Quality evaluation index of shipwreck image

從表1 中的數(shù)據(jù)可以明顯看出，采用小波變換（WT）、脊波變換（RT）和輪廓波變換（CT）去噪后，圖像的SNR 值和PSNR 值有一定程度的提升，但SSIM值效果不理想，圖像失真較為嚴(yán)重.采用剪切波變換（ST）去噪后，SNR 值、PSNR 值和SSIM 值分別提升了232.87%、59.85%和60.73%，MSE 值下降了86.46%.采用非下采樣剪切波變換（NSST）去噪后，SNR 值、PSNR 值和SSIM 值分別提升了253.05%、65.04%和65.39%，MSE 值下降了88.40%，去噪效果好于剪切波變換.運用密度聚類、灰度變換和非下采樣剪切波變換相結(jié)合的方法處理后，圖像的SNR 值從3.702 8 dB 提升至13.855 7 dB，提升了274.20%；PSNR 值從14.406 8 dB 提升至24.559 7 dB，提升了70.47%；SSIM 值從0.328 2 提升至0.565 2，提升了72.21%；MSE 值從0.036 2 下降至0.003 5，下降了90.33%.本文方法的SNR 值、PSNR 值和SSIM 值提升的幅度以及MSE值下降的幅度均大于小波變換（WT）、脊波變換（RT）、輪廓波變換（CT）、剪切波變換（ST）和非下采樣剪切波變換（NSST）去噪法，從客觀數(shù)據(jù)表明，本文方法在上述6 種變換域去噪法中去噪效果最好，失真程度最小，邊緣保持能力最強(qiáng).

實驗2為進(jìn)一步驗證本文方法在不同噪聲密度條件下的客觀性和優(yōu)越性，利用公共圖像數(shù)據(jù)集進(jìn)行實驗，并與現(xiàn)有文獻(xiàn)方法進(jìn)行對比.將圖9 中的公共圖像數(shù)據(jù)集添加方差分別為0.3、0.4、0.5 和0.6的斑點噪聲，運用文獻(xiàn)［25］的BM3D 算法、文獻(xiàn)［26］的PGPD 算法、文獻(xiàn)［27］的EPLL-GMM 算法、文獻(xiàn)［28］的ACVA 算法和本文所提方法進(jìn)行去噪處理，得到去噪后圖像的MSE、SNR、PSNR 和SSIM 的平均值如表2所示.

圖9 公共圖像數(shù)據(jù)集Fig.9 Public image dataset

表2 公共圖像數(shù)據(jù)集質(zhì)量評價指標(biāo)Tab.2 Quality evaluation index of public image dataset

從表2 中的數(shù)據(jù)可知，本文所提方法在噪聲方差分別為0.3、0.4、0.5和0.6的情況下，MSE 值分別下降了89.76%、91.21%、91.40% 和91.64%；SNR 值、PSNR 值和SSIM 值分別提升了164.34%、197.74%、223.60%、247.41%，86.42%、97.90%、103.25%、110.20%和293.32%、345.45%、387.73%、422.05%.各項指標(biāo)的變化幅度隨噪聲方差的增大而逐漸提升，并且其值均優(yōu)于文獻(xiàn)［25-28］的去噪方法，表明本文方法在去除聲吶圖像斑點噪聲中擁有明顯優(yōu)勢.通過折線圖進(jìn)一步展示各項評價指標(biāo)的變化趨勢，趨勢圖如圖10所示.

圖10 圖像質(zhì)量評價指標(biāo)趨勢圖Fig.10 Trend charts of image quality evaluation indexes

圖10 中，隨著噪聲方差的增大，BM3D、PGPD、EPLL-GMM、ACVA 和本文所提方法的MSE 值呈上升趨勢，SNR 值、PSNR 值和SSIM 值呈下降趨勢.其中，BM3D 的MSE 值在5 種去噪方法中保持最大，SNR值和PSNR 值保持最小，表明其去噪效果不太理想；ACVA 的SNR 值、PSNR 值 和SSIM 值均大于PGPD 和EPLL-GMM 算法，表明其去噪效果好于PGPD 和EPLL-GMM 算法；本文去噪法的各項指標(biāo)變化平穩(wěn)，始終優(yōu)于其他4 種去噪方法，并且其優(yōu)勢隨噪聲方差的增大而更加明顯，證明了本文方法更適用于具有高密度噪聲的聲吶圖像去噪.不同方法的去噪效果如圖11所示.

圖11 不同方法的去噪效果Fig.11 Denoising effect of different methods

在圖11 中，含噪圖像分別通過BM3D、PGPD、EPLL-GMM、ACVA 和本文方法進(jìn)行去噪處理后，圖像質(zhì)量得到了不同程度的提升.其中BM3D 算法去除了圖像中的部分斑點噪聲，但效果還不夠理想；PGPD 和EPLL-GMM 算法噪聲去除明顯，但細(xì)節(jié)信息丟失較為嚴(yán)重，影響了圖像的質(zhì)量；ACVA 算法的邊緣保持能力得到較大的提升，但與本文方法相比還存在一定差距.

除了MSE、SNR、PSNR 和SSIM 等常用指標(biāo)，時間復(fù)雜度也是評價圖像去噪算法的重要標(biāo)準(zhǔn).運用搭載Intel Core i5-6200 處理器和4 GB 運行內(nèi)存的計算機(jī)在MATLAB R2019b 平臺上進(jìn)行多次試驗，得到不同算法的平均運行時間如表3所示.

表3 不同算法的平均運行時間Tab.3 Average operation time of different algorithms

由表3可知，BM3D的平均運行時間為0.853 8 s，時間復(fù)雜度低，但卻以犧牲去噪性能為代價；本文方法平均運行時間為2.594 1 s，明顯低于PGPD、EPLLGMM 和ACVA 去噪算法，驗證了本文方法較強(qiáng)的去噪實時性.

4 結(jié)論

本文對聲吶圖像去噪方法進(jìn)行研究，提出一種基于密度聚類與灰度變換的非下采樣剪切波域圖像去噪方法，該研究可總結(jié)為：

1）利用非下采樣剪切波變換的局部化特性、方向敏感性、最優(yōu)稀疏表示特性和平移不變性等特點降低去噪后圖像的失真程度，緩解偽吉布斯現(xiàn)象，降低了時間復(fù)雜度；

2）對聲吶圖像的噪聲特性進(jìn)行分析，選擇適合其分布的DBSCAN 聚類算法，有效分離圖像中的噪聲干擾，增強(qiáng)圖像的邊緣保持能力；

3）利用灰度變換中的分段函數(shù)增強(qiáng)圖像的對比度，提升了視覺效果；

4）仿真試驗表明，本文方法在均方誤差、峰值信噪比和結(jié)構(gòu)相似度等指標(biāo)上優(yōu)于傳統(tǒng)變換域和文獻(xiàn)［25-28］中的去噪方法，并且在高密度噪聲條件下其優(yōu)勢更為明顯，適用于噪聲干擾嚴(yán)重的聲吶圖像去噪.