楊 培，高雷阜，訾玲玲

1.遼寧工程技術大學 理學院，遼寧 阜新 123000

2.遼寧工程技術大學 電子與信息工程學院，遼寧 葫蘆島 125105

1 引言

圖像融合是圖像處理的關鍵技術之一，在軍事、遙感、機器人等領域廣泛應用，其目的是提高信息的利用率，最大限度提取同一目標的多源圖像數(shù)據(jù)綜合成高質量圖像。為此，學者們相繼提出許多圖像融合算法[1-3]，其中基于多尺度變換的方法，如小波變換[4]、輪廓波變換[5]、非下采樣輪廓波變換（nonsubsampled contourlet transform，NSCT）[6]等，以其多尺度、多方向、各向異性以及平移不變等優(yōu)點被廣泛應用。但經(jīng)這些多尺度變換后的低頻圖像系數(shù)不具備稀疏性，對圖像的特征不能準確逼近。而稀疏表示（sparse representation，SR）是通過對源圖像的分解建立細節(jié)特征和稀疏系數(shù)的映射關系，提取出圖像更深層次的結構信息。因而，學者們將稀疏表示理論引入圖像融合領域來解決多尺度變換的特征逼近問題[7-8]。Hu等[9]提出了基于SR和IHS（intensity,hue,saturation）變換的遙感圖像融合，先對圖像進行IHS變換，再用離散余弦變換（discrete cosine transform，DCT）字典進行稀疏表示，然后用絕對值最大融合規(guī)則對稀疏系數(shù)進行融合，最后通過重構和IHS逆變換得到融合圖像。雖然達到了一定的融合效果，但其中的DCT字典進行稀疏表示時丟失了邊界信息，絕對值取大的融合規(guī)則也易丟失圖像間的結構信息。為此，嚴春滿等[10]采用滑窗模型和K-SVD（K-singular value decomposition）字典學習算法對多聚焦圖像進行融合，解決了稀疏表示的邊界結構問題，但耗時較長。劉慧等[11]用NSCT變換和IHS變換進行結合，既克服了頻譜失真問題，又保證了多方向性和平移不變性，但算法復雜度很高。趙春暉等[12]提出了一種快速的基于SR和NSCT變換的圖像融合算法，稀疏表示采用4方向法（four-direction sparse representation，F(xiàn)DSR）提升融合效率，但在融合規(guī)則不恰當時易產(chǎn)生斜紋效果。

針對以上問題，本文提出了螺旋結構和梯度分析的圖像融合算法。將待融合圖像進行NSCT變換，對低頻信息，利用螺旋處理前后的“2+2”方向結構圖像數(shù)據(jù)進行字典訓練和稀疏表示，既避免了傳統(tǒng)滑動窗口帶來的冗余計算量，又消除了FDSR的斜紋效果，再利用空間頻率信息進行系數(shù)融合；對高頻信息，先對其進行梯度分析，保留圖像的結構信息，再采用最大模值優(yōu)先進行融合；最后通過實驗驗證本文方法的有效性。

2 基礎理論

2.1 NSCT變換理論

NSCT變換是由非下采樣金字塔濾波器（nonsubsampled pyramid filter bank，NSPFB）和非下采樣方向濾波器組（nonsubsampled directional filter bank，NSDFB）組成。NSCT變換如圖1所示[6]。其中，NSPFB保證了多尺度特性，其特性是由移位不變?yōu)V波結構獲得的，通過雙通道非采樣二維濾波器組實現(xiàn)了拉普拉斯算子的金字塔子帶分解。NSDFB給出了方向性，將臨界采樣雙通道扇形濾波器組和重采樣相結合，構建樹形結構的濾波器組，將二維頻率平面分割成方向楔形，利用NSDFB得到位移不變的方向展開。

NSCT變換能夠很好地分解出圖像中的細節(jié)和近似光譜信息，而且在圖像的分解和重構過程中，避免了對圖像的上下重采樣，使得各個子帶圖像與源圖像的尺寸大小保持一致，解決了頻譜混疊效應，使所表達的圖像更精細。

Fig.1 Nonsubsampled contourlet transform圖1 NSCT變換

2.2 稀疏表示和字典學習

稀疏表示是用字典中少量原子的線性組合來描述信號。設信號x∈RN，字典D∈RN×L，稀疏系數(shù)α∈RL，則稀疏表示問題可以描述為：

其中，||?||0表示l0范數(shù)，表示非零元素的個數(shù)。這是一個NP難問題。Donoho和Candes等[13-14]已經(jīng)證明在一定的條件下，l0范數(shù)可以轉換為l1范數(shù)優(yōu)化問題。進而，式（1）采用正則化方法可得到如下表達[15]：

其中，λ為正則化參數(shù)，用來均衡稀疏性和精度。這樣的正則化稀疏模型有很多求解算法，最小角回歸（least angle regression，LAR）算法是其中一種高效的解法[16]。

在圖像的稀疏表示中，字典D的選取會影響到分解的稀疏度和重構的精度。構造字典主要有兩種方法：解析方法和學習方法。解析方法由數(shù)學模型計算得到字典，獲取簡單且普適性高，但稀疏度和誤差也高，如DCT字典等。學習方法通過樣本數(shù)據(jù)進行字典訓練，誤差和稀疏度相對較小，但訓練成本較高且普適性稍差，如K-SVD算法、ODL（online dictionary learning）[17]算法。K-SVD實時性較差，字典迭代更新時ODL的效率更高，因而本文采用學習方法在圖像融合過程中進行字典訓練，ODL算法進行字典學習，LAR算法計算稀疏編碼，兩者交替迭代，快速準確提取圖像特征。

稀疏表示是在單尺度上對圖像信息進行分析，與采用多尺度分析的NSCT變換進行結合，能滿足多尺度上提取圖像重要特征的要求。通過NSCT變換獲得的圖像低頻信息，直接對其進行融合處理會導致無法獲得圖像中顯著特征信息，并且變換中的多尺度分析需對圖像進行分解和合成，易導致樣本信息丟失，而稀疏表示中字典學習算法能對圖像中的重要特征進行特征提取，恰好可以克服多尺度分析中的不足，并且字典原子是通過待融合圖像信息訓練而來，可以獲得圖像更準確的結構描述。故而本文在NSCT變換的多尺度基礎上，對圖像低頻子圖進行稀疏逼近表示再融合處理，對已經(jīng)稀疏化的高頻子圖，通過梯度分析讓邊緣特征更加凸顯后再融合，使最終獲得的融合圖像有更優(yōu)的效果。

3 螺旋結構及梯度分析的融合算法

本文算法由兩部分組成：第一部分是NSCT變換后低頻子圖中的螺旋結構設計；第二部分是高頻子圖的融合規(guī)則。在低頻子圖中采用螺旋結構，即是“2+2”方向稀疏表示進行字典學習、稀疏融合。高頻子圖中進行梯度分析后再進行融合。整體過程如圖2所示。

3.1 低頻融合過程

3.1.1 螺旋結構設計

稀疏表示融合規(guī)則中，圖像經(jīng)常被劃分為大小相同的塊，采用滑動窗口的形式對每個圖像塊進行稀疏表示和融合，滑動窗口的步長越小，融合效果越好。這種處理方式中，若融合規(guī)則不匹配，易出現(xiàn)空窗現(xiàn)象。此外，步長越小則計算量越大，在n×n大小的窗口滑動過程中，一個像素信息在步長為1時會被重復使用n2次，這極大降低了融合效率。文獻[12]中放棄滑窗模型，用FDSR進行圖像稀疏表示，重復計算次數(shù)從n2次降低到了4次且避免了空窗現(xiàn)象。雖然提供了更多的方向信息，但實驗中發(fā)現(xiàn)這種處理融合圖像易出現(xiàn)斜紋效果。對于二維圖像來說，空間基底應是兩個線性無關的向量。正交的向量組構成的空間更是適應人類視覺的感受。為了解決上述問題并且考慮空間圖像的基底及其旋轉前后的正交設計要求，本文提出了一種螺旋結構設計方法，也稱之為螺旋處理前后的“2+2”方向稀疏表示（double_two direction sparse representation，DtDSR）。對源圖像在螺旋處理前進行兩個方向（45°方向和135°方向）稀疏表示，螺旋處理后又進行兩個方向（水平方向和豎直方向）稀疏表示，具體過程如圖3所示，即：

Fig.2 Flow diagram of proposed algorithm圖2 本文算法流程框架圖

其中，A0、B0分別為圖像A、B的低頻子圖。

對于圖像中任一像素點（邊緣點除外）來說，周邊有8個像素是其最近鄰相關點，螺旋分析前的兩個方向，保證了其對角線方向的4個像素依舊與原像素點最近鄰，螺旋分析后所采用的兩個方向，也較大限度拓撲同構了原像素點周邊的結構性。此外，稀疏表示是對零范數(shù)NP難問題的凸優(yōu)化處理，以便于利用相對成熟的凸優(yōu)化理論和方法進行求解，故而采用稀疏字典對每個方向上的列或行像素進行稀疏表示以得到最優(yōu)解，再進行系數(shù)融合。又因4個方向所獲得的樣本在處理過程中對任一像素均只使用一次，所占權重相同，因而在最終融合時，對所獲得的4個不同方向上的融合圖像進行平均得到低頻融合結果。

在DtDSR結構模型中，每個像素信息僅被重復使用4次，相比于滑動窗口模型的每個像素被重復使用n2次（n為滑動窗口的長度）來說，所消耗的時間會大幅度下降，這也極大提高了融合算法的效率。

3.1.2 低頻融合過程具體步驟

步驟1將低頻子圖A0、B0分別進行邊緣擴展，避免后期稀疏表示出現(xiàn)邊緣效應并進行45°、135°方向處理，各自得到訓練樣本A1、A2和B1、B2；將A0、B0分別進行螺旋結構處理后再進行水平、豎直方向處理，各自得到兩個訓練樣本A3、A4以及B3、B4；具體處理過程如3.1.1小節(jié)所述。

步驟2對兩組樣本A1、A2、A3、A4和B1、B2、B3、B4根據(jù)式（2）分別進行ODL字典學習、LARs稀疏編碼，得最優(yōu)字典和稀疏系數(shù)，即：

步驟3對稀疏系數(shù)，采用文獻[7]的系數(shù)融合方法，分別計算其空間頻率，取值較大的作為融合系數(shù)，并同時選定對應的字典作為融合后的字典

步驟4重構圖像，即：

其中，F(xiàn)表示最終的低頻融合圖像。

3.2 高頻融合過程

Fig.3 Schematic diagram of DtDSR structural model圖3 DtDSR結構模型示意圖

NSCT分解所得的高頻子圖含有圖像的邊緣和細節(jié)等突變特性。高頻系數(shù)絕對值的大小反映這些特征變化的劇烈程度。對于圖像的平滑部分，高頻子圖系數(shù)絕對值比較??；對于圖像的高頻邊緣或者細節(jié)突變部分，高頻子圖系數(shù)的絕對值一般比較大。通常情況下，對于高頻子圖的融合常采用“模值取大”這一策略，但對于某些待融合圖像來說，在同一部分都有突變，并且高頻子圖系數(shù)值都比較大，此時，選取到的高頻子圖系數(shù)容易使融合圖像出現(xiàn)模糊現(xiàn)象。如圖4所示，兩張待融合圖片分別為計算機斷層成像得到的彩色圖片A和核磁共振得到的灰度圖像B，在有些邊緣處如圖中紅框區(qū)域，視覺觀察A較B模糊，但實驗中，經(jīng)過NSCT變換后得到的高頻子圖中，A系數(shù)的模值大于B的值，導致融合時會選取較模糊的A。為了減少這種現(xiàn)象的出現(xiàn)并且保持融合效率，提出一種基于梯度分析（gradient analysis，GA）的高頻融合方法。

在圖像處理中，圖像的邊緣可通過計算梯度求得，所得結果反映像素點灰度的最大變化率及其方向，能準確判斷出圖像邊緣。若像素點處于圖像邊緣，則具有較大梯度值；若處于非邊緣區(qū)域，梯度值很小或者接近零。此外，梯度運算也是圖像邊緣特征提取的高效方法之一，根據(jù)梯度的這些特性，將高頻子圖的邊緣梯度表示引入“模值最大”這一融合策略中。如上例中A、B的高頻子圖，由于A的模糊性導致圖中紅框內邊緣變化不明顯，在同一邊緣處經(jīng)過梯度分析后得到新的系數(shù)后，A的梯度分析系數(shù)小于邊緣比較清晰的B的值，融合時會選擇B，這樣會達到更好的融合效果。高頻子圖融合規(guī)則具體過程如下所示。

步驟1對A的高頻子圖，求其水平和豎直方向上的梯度表示矩陣，即：

Fig.4 Part of obtained graph data before and after gradient analysis圖4 梯度分析前后部分數(shù)據(jù)圖

步驟3根據(jù)圖像A、B各尺度、各方向的梯度系數(shù)矩陣決定融合圖像F的高頻選取策略，如下式所示：

37 Attention of community population in Shanghai to chronic kidney disease and related influencing factors

3.3 融合過程算法步驟

輸入：圖像A、B。

輸出：融合圖像IF。

（1）對A、B進行NSCT分解，分解級數(shù)根據(jù)圖像大小采用自適應方式進行選取，得到低頻子圖A0、B0和高頻子圖

（2）低頻處理過程。由低頻子圖A0、B0，根據(jù)式（3）得到低頻處理的訓練樣本，然后用式（4）得到各自的字典和稀疏系數(shù)。采用空間頻率取大的方式得到融合系數(shù)和融合字典。再利用式（5）和式（6）得到最終低頻融合圖像F。

（4）將獲得的F和進行NSCT反變換，獲得最終的融合圖像IF。

4 實驗及結果分析

本文算法定義為DtDSR-GA，為了驗證融合算法的性能，采用另外五種融合算法進行對比實驗?；赟R的圖像融合方法[18]（SR）、基于小波變換和稀疏表示的圖像融合算法[7]（DWT-SR）、基于非下采樣輪廓波變換和稀疏表示的融合算法[19]（NSCT-SR）這三種算法中融合規(guī)則選取按照相應文獻中的描述設定。此外，在NSCT變換的基礎上，低頻采用不同方向結構（FDSR、DtDSR）稀疏表示，融合規(guī)則同為低頻分量稀疏系數(shù)空間頻率取大，高頻分量“模值取大”，算法分別記為NSCT-FDSR、NSCT-DtDSR。

為了能更好地凸顯本文算法的效果，實驗在相同基礎條件下進行比較，五種對比算法和本文算法在稀疏表示時均采用ODL字典學習和LARs稀疏編碼[20]，正則化參數(shù)為0.15，迭代次數(shù)為10，滑動窗口大小為8×8，步長為1。DWT變換時小波基使用“db4”，分解層數(shù)設置為4層。NSCT變換時金字塔濾波器組選用“maxflat”，方向濾波器組采用“dmaxflat7”，分解尺度統(tǒng)一采用自適應于圖像尺寸的方式，分解方向設置為8。本文算法利用Matlab R2014a編程實現(xiàn)，運行環(huán)境Intel?Core?i3-6100U CPU，4 GB內存，操作系統(tǒng)Windows 7。

4.1 評價指標

融合圖像的質量評價是評估融合方法性能的重要手段，實驗在主觀評價的基礎上，用客觀評價作為輔助手段來驗證算法的有效性。本文選取以下5個指標作為客觀評價的標準[21]。在以下評價指標的公式中，F(xiàn)表示融合圖像，R表示參考圖像，M、N是圖像F的大小，F(xiàn)(i,j)為F在(i,j)點的像素值。

（2）標準差（standard deviation，SD）

SD衡量圖像灰度的離散程度。設是圖像的平均值，則SD的計算公式為：

（3）空間頻率（spatial frequency，SF）

SF表示圖像中的空間像素總體活躍度。設fx、fy分別是x方向的行頻率和y方向的列頻率，則SF的計算公式為：

（4）平均梯度（average gradient，AG）

AG體現(xiàn)圖像中的微小細節(jié)反差和紋理變化特征，其計算公式為：

（5）平均結構相似度（average structural similarity，ASSIM）

結構相似度SSIM[22]是圖像結構保真度的度量，其計算公式為：

其中，l(R,F)、c(R,F)、s(R,F)分別是亮度、對比度、結構函數(shù)。uR和uF為平均值，σR和σF為標準差，σRF為協(xié)方差，C1、C2和C3是為了避免分母接近零而產(chǎn)生不穩(wěn)定所添加的非常小的正常數(shù)。結構相似度越大，融合效果越好。為了減少評價誤差，采取平均結構相似度ASSIM進行評價，其定義如下：

以上五種評價指標都是越大越好，對比實驗表中最優(yōu)指標由粗體標示。

4.2 實驗待融合圖像集

實驗分別對灰度圖像、彩色圖像進行仿真實驗，所采用圖像均經(jīng)過幾何校正、圖像配準等預處理步驟。實驗數(shù)據(jù)包括“Disk”（256×256）和“Pepsi”（512×512）兩組灰度圖像集，以及“Street”（600×475）、“City”（256×256）、“Medical”（256×256）、“IR-visible”（640×480）四組圖像集，如圖5所示。

4.3 灰度圖像融合實驗

圖6給出了灰度圖像主觀視覺融合效果圖，所采用融合方法從左至右依次是SR、DWT-SR、NSCT-SR、NSCT-FSR、NSCT-DtDSR以及本文方法DtDSR-GA。觀察圖6可以看出，SR方法融合效果有模糊現(xiàn)象，DWT-SR方法所獲取融合圖像的清晰度稍遜NSCT系列方法和本文方法。尤其在Disk效果圖中，本文方法明顯要優(yōu)于其他融合方法。

將Pepsi圖像集及融合效果圖中箭頭和“P”中間部分放大，如圖7所示，融合算法次序同圖6一致。Pepsi圖像集中白色偏斜紋，本文算法很好地在融合圖像中繼承下來，且圖中箭頭部分的字跡比較清晰。SR方法由于其融合規(guī)則的不足，融合圖像清晰程度較低。DWT-SR融合方法因采用下采樣操作，出現(xiàn)頻譜混疊現(xiàn)象，在融合結果中表現(xiàn)出模糊鋸齒狀。NSCT-SR方法雖使得圖像邊緣變得清晰，但較本文算法依然有差距。NSCT-FSR、NSCT-DtDSR在圖7中與本文算法效果沒有過于明顯的區(qū)別。

表1是灰度圖像集融合結果的五種客觀評價指標對比。由表1中可以看出，用本文提出的方法DtDSR-GA得到的融合圖像，五項指標較其他算法都有提升。其信息熵IE、空間頻率SF的提高，表明融合圖像保留了源圖像的大部分信息，輪廓清晰；標準差SD的提高，表明融合圖像灰度變化劇烈，對比度較好；平均梯度AG的提高，表明了融合結果細節(jié)清晰程度。本文算法與NSCT-FSR、NSCT-DtDSR在圖7中相比較沒明顯區(qū)別，但在表1中的數(shù)據(jù)顯示其依然有較高提升。

Fig.5 Image data sets圖5 圖像數(shù)據(jù)集

Fig.6 Fusion result of grayscale images圖6 灰度圖像融合效果圖

Fig.7 Local enlarged images of fusion effect in“Pepsi”圖7 “Pepsi”融合效果局部放大圖

Table 1 Comparison of evaluation indices of fusion results of grayscale image sets表1 灰度圖像集融合結果評價指標對比

4.4 彩色圖像融合實驗

本文融合算法對于彩色圖像或彩色與灰度圖像之間的融合也同樣適用，其流程如圖8所示。

Fig.8 Flow chart for color image fusion algorithm圖8 彩色圖像融合流程圖

實驗的主觀效果如圖9所示，融合算法次序同圖6一致?！癝treet”是多聚焦彩色圖像，本文算法的融合效果明顯優(yōu)于其他算法。“City”是遙感圖像，融合結果中采用NSCT分解后的視覺效果較優(yōu)，本文算法在紋理結構上比其他算法更為清晰?！癕edical”中彩色圖像是顯示人體組織血流和代謝功能的計算機斷層成像圖像，灰度圖像是對軟組織和血管清晰的分辨率高的核磁共振圖像。融合效果圖中，SR與NSCTSR算法色彩效果近似，但與源圖像有較大差距。DWT-SR色彩效果同源圖像近似，但結構信息模糊。本文算法融合效果不但在色彩上接近源圖像，在紋理結構上也達到了較好的效果。“IR-visible”是可見光和紅外圖像，融合的色彩效果與源可見光圖主觀比較無色差。對于紅外圖像的煙霧彈亮點以及煙霧掩蓋下的結構紋理部分，本文算法融合后效果更好。且通過表2客觀指標對比，本文算法不但主觀上效果良好，且客觀指標較其他比較算法也都有所提高。

4.5 算法時效實驗

Fig.9 Fusion result of color images圖9 彩色圖像融合效果圖

Table 2 Comparison of evaluation indices of fusion results of color image sets表2 彩色圖像集融合結果評價指標對比

Table 3 Comparison of running time of different fusion methods for different image sets表3 不同圖像集的不同融合方法運行時間比較 s

圖像融合中，時效性也是評價算法的重要指標。本文算法在時間上的優(yōu)勢體現(xiàn)在低頻融合時的DtDSR結構模型上，該模型中對圖像的任一像素總共只使用了4次，而滑窗模型在步長為1，窗口大小為n×n時，對每個像素使用了n2次，樣本數(shù)據(jù)量的減少使得本文算法復雜度在數(shù)量級上有較大降低，在時間上占有一定的優(yōu)勢。表3給出了實驗中不同算法的融合時間。SR、DWT-SR、NSCT-SR方法采用滑窗模型，在步長為1，窗口大小為8時每個元素重復使用了64次，導致融合時間比較長。DWT-SR由于圖像下采樣，融合時間相對另兩種方法有所降低。而NSCT-FSR、NSCT-DtDSR、DtDSR-GA采用方向模型，每個元素只用4次，提高了融合效率。NSCTDtDSR方法因增加了螺旋處理過程，運行時間比NSCT-FSR稍長，DtDSR-GA與NSCT-DtDSR相比，在高頻處理中由于增加了梯度運算過程使得時間也稍增。與傳統(tǒng)的滑窗模型的稀疏融合算法相比較，本文提出的DtDSR-GA算法數(shù)十倍提高了融合速度，雖在高頻子圖中因梯度分析使得時間有所增加，但從整體融合時間來說，增加時間與總時間的比值不超過3%。結合前文的主客觀融合質量分析，DtDSRGA不但融合效率高，且融合質量也有所提高，證明了其有效性。

5 結束語

本文提出一種基于NSCT變換的圖像融合方法，在低頻子圖融合中，提出一種螺旋處理的“2+2”方向結構模型進行字典學習和稀疏分解。在高頻子圖融合過程中，先對高頻分量進行梯度分析處理后再采用“模值取大”的融合規(guī)則。通過實驗證明，采用本文方法的融合結果在主觀和客觀上都具有一定的優(yōu)勢。本文的研究工作還可以有進一步的改進和應用，如，融合規(guī)則中可以采用多特征處理融合，采用本文稀疏處理方法嘗試對彩色圖像和RGB-D圖像進行融合等。