李 隆，傅依柳，柴昱洲，陳霄鵬，高當(dāng)麗，屈子杰

引言

圖像融合是將兩幅或兩幅以上的源圖像組合成一幅包含源圖像顯著信息的復(fù)合圖像。一般來說圖像融合的目的就是提取源圖像中的有用信息，并將其結(jié)合、保留到融合圖像中。對于不同的應(yīng)用場合、不同的圖像源，其融合策略也不盡相同。20世紀(jì)90年代開始，隨著小波理論的產(chǎn)生和不斷完善，這種具有良好時(shí)頻定位特性的方法在圖像融合領(lǐng)域得到了極為廣泛的應(yīng)用，從而開啟了多尺度分解應(yīng)用于圖像融合的大門。然而小波變換只能捕捉水平、垂直和對角3個(gè)方向的信息，具有較差的方向性。針對小波變換在圖像處理中的缺陷，Do和 Vetteri［1］提出了一種“真正”的二維圖像稀疏表達(dá)方法——contourlet變換。它不僅具有多尺度和良好的時(shí)頻局部特性，還具有多方向性和各向異性。但是，在對圖像進(jìn)行contourlet變換的過程中，需要對圖像進(jìn)行下采樣操作，這使得contourlet變換不具備平移不變性，導(dǎo)致在圖像處理中產(chǎn)生偽吉布斯現(xiàn)象。為此 A.L.Cunha等［2］又提出了具有平移不變性的contourlet變換——非采 樣 contourlet變 換 （NSCT），但 是NSCT計(jì)算復(fù)雜度較高且運(yùn)行時(shí)間較長。近年來，K.Guo等［3］提出了一種更優(yōu)的多尺度幾何分析工具——shearlet變換（ST），這種變換結(jié)構(gòu)簡單并且易于實(shí)現(xiàn)，能對圖像中的邊界輪廓產(chǎn)生最優(yōu)逼近。非下采樣剪切波變換［4］（NSST）作為ST的改進(jìn)型模型，具有優(yōu)越的圖像處理性能，避免了圖像采樣中的細(xì)節(jié)丟失，并且具有平移不變性，但在圖像融合領(lǐng)域的應(yīng)用仍處于初步探索階段。

對于分解后子帶圖像的融合，目前已有的融合策略中，低頻子帶多采用加權(quán)平均的方法，融合系數(shù)是各源圖像的局部或全局顯著特征（如均值，標(biāo)準(zhǔn)差）結(jié)合不同權(quán)重所得；高頻子帶多采用區(qū)域窗口的顯著特征選取系數(shù)或加權(quán)系數(shù)的方法。上述方法在融合過程中易受噪聲干擾，導(dǎo)致融合結(jié)果引入虛假信息。為了改善上述缺陷，文獻(xiàn)［5］提出了一種基于改進(jìn)梯度投影非負(fù)矩陣分解（NMF）融合策略。非負(fù)矩陣分解是一種新的矩陣分析方法，它能夠更好地提取和描述圖像局部的特征信息，可實(shí)現(xiàn)真正的局部特性提取，從而得到較好的融合效果。因此本文在低頻系數(shù)融合策略方面采用NMF分解，對于高頻子帶的融合，通過計(jì)算區(qū)域窗口的高斯加權(quán)區(qū)域能量及區(qū)域方差，然后分別歸一化結(jié)合給定的加權(quán)系數(shù)，得到融合后的高頻子帶系數(shù)。

1 基本原理

1.1 NSST相關(guān)理論

根據(jù)合成小波的相關(guān)理論，剪切波變換是通過仿射系統(tǒng)將幾何和多尺度相結(jié)合起來，當(dāng)維數(shù)n＝2時(shí)，具有合成膨脹的仿射系統(tǒng)表達(dá)式為

式中：ψ∈L2（R2）；L2（R2）表示平方可積空間；A、B均為2×2的可逆矩陣，且｜det（B）｜＝1。若ΨAB（ψ）中元素滿足Parseval框架時(shí)，則ΨAB（ψ）的元素稱為合成小波，A稱為各向異性膨脹矩陣，Aj與尺度變換相關(guān)聯(lián)；B為剪切矩陣，Bl與幾何變換相關(guān)聯(lián)，如剪切和旋轉(zhuǎn)變換。由（1）式可以看出，合成小波具有多分辨率的特性。

K.Guo等人提出shearlet變換系統(tǒng)如下。通過對ψ∈L2（R2）進(jìn)行膨脹剪切和平移變換形成剪切波：

shearlet能夠充分表征物體細(xì)節(jié)特征，是由于其支撐基具有一定斜率且成梯形分布，根據(jù)切波理論，其支撐基是一對尺寸近似為22j×2j，斜率2－jL的梯形支撐區(qū)間（j為分解尺度，L為方向因數(shù)），如圖1所示。

圖1 NSST多尺度多方向分解過程Fig.1 NSST decomposition in multi－scale and multi－direction

同NSCT的結(jié)構(gòu)類似，NSST的離散化過程同樣分為多尺度分解和多方向分解2部分，其中多尺度分解仍為NSCT中的非下采樣金字塔分解方式，源圖像經(jīng)R級非下采樣濾波器分解可產(chǎn)生1個(gè)低頻子帶圖像和R個(gè)高頻子帶圖像，分解后的圖像大小與源圖像大小相同；其次NSST的多方向分解是通過改進(jìn)的剪切波濾波器實(shí)現(xiàn)，即把標(biāo)準(zhǔn)的剪切濾波器從偽極化網(wǎng)絡(luò)系統(tǒng)映射回笛卡爾坐標(biāo)系統(tǒng)，經(jīng)過L級方向分解可得到n＝2L＋2個(gè)方向子帶Jjn。如圖2所示。

圖2 NSST梯形支撐區(qū)間Fig.2 Range of NSST trapezoidal support

1.2 改進(jìn)梯度投影非負(fù)矩陣分解

非負(fù)矩陣分解［6］（NMF）是一種新的矩陣分析方法，它可以描述為對于一個(gè)非負(fù)矩陣V，求非負(fù)的矩陣W∈Rn×r＋和非負(fù)的矩陣H∈Rr×m＋，使V≈WH，分解后的矩陣僅包含非負(fù)元素，只允許獲得的數(shù)據(jù)是原始數(shù)據(jù)的加性組合，從而使原矩陣中列向量可解。

分解前后可以理解為：原始矩陣V中的列向量是對左矩陣W 中所有行向量（稱為基向量）的加權(quán)和，而權(quán)重系數(shù)就是右矩陣中對應(yīng)列向量中的元素。這相當(dāng)于用相對少的基向量來表示大量的數(shù)據(jù)向量，算法所得到的非負(fù)基向量組W 具有一定的線性無關(guān)性和稀疏性，從而使其對原始數(shù)據(jù)的特征及結(jié)構(gòu)具有相當(dāng)?shù)谋磉_(dá)能力。因此，只要適當(dāng)?shù)剡x取基向量的維數(shù)，就可以獲取原始數(shù)據(jù)的特征。

NMF的求解問題實(shí)際上是一個(gè)最優(yōu)化問題，利用乘性迭代的思路求解W 和H。為了降低NMF算法的時(shí)間復(fù)雜度，這里用改進(jìn)梯度投影的方法［5］。NMF的求解可改寫為2個(gè)目標(biāo)函數(shù)，即：

根據(jù)（3）～（4）式，對于當(dāng)前結(jié)果（ˉW，ˉH），同時(shí)使用梯度投影算子：

其中P（x）＝max（0，x），α為步長。更新α的算法如下：

步驟1）輸入非負(fù)矩陣V，初始化W1和H1；

步驟2）給定0＜β＜1，0＜σ＜1，設(shè)α0＝1，β＝0.1，β＝0.1，σ＝0.01；

步驟3）對于k＝1，2，…，αk←αk－1，若αk滿足：

其中αk＝βtk，tk是非負(fù)整數(shù)t的第一項(xiàng)，且

則反復(fù)執(zhí)行αk←αk／β，直到不滿足（5）式或 H（αk／β）＝H（αk）時(shí)停止。若αk不滿足（5）式，則執(zhí)行αk←αk·β直到滿足（5）式。然后設(shè) Hk＋1＝P［Hkαk▽f（Hk）］；

步驟4）若Hk＋1滿足（5）式，Wk＋1的計(jì)算與上面過程類似則輸出W 和H，否則重復(fù)步驟2）～3）。

從上述過程可看出：ˉf（?W，?H）＝1／2·‖V－?W?H‖2F的時(shí)間復(fù)雜度僅為O（nmr），從而使整體時(shí)間復(fù)雜度為O（tnmr），t是（5）式中每次循環(huán)的平均值。

2 圖像融合策略

圖像的多尺度融合主要包括多尺度分解、高頻低頻子帶融合以及圖像重構(gòu)3個(gè)過程。本文實(shí)驗(yàn)數(shù)據(jù)為兩幅圖像的融合，設(shè)源圖像A、B經(jīng)過NSST分解后得到不同尺度的方向子帶系數(shù)｛LAJ，HAj，k｝、｛LBJ，HBj，k｝。其中J 表示分解層數(shù)，LXJ（X＝A，B）表示圖像X的低頻子帶系數(shù)，HXj，k表示圖像X 的第j層第k個(gè)高頻子帶系數(shù)［8］。

2.1 低頻子帶融合策略

圖像融合中，參與融合的觀測圖像本質(zhì)上就是真實(shí)圖像通過不同傳感器并且經(jīng)過加噪得到的。低頻子帶是在低分辨率下對原始圖像的近似，包含了圖像的大量基本信息，為了最大程度地綜合圖像的特征信息，這里采用1.2節(jié)提到的算法。

設(shè)觀測圖像為f1，f2，…，fk（k≥2），圖像像素大小均為m×n，將觀測圖像的逐個(gè)元素按照優(yōu)先的方式儲存到一個(gè)列向量中，得到k個(gè)列向量v1，v2，…，vk，其中將這k個(gè)列向量的數(shù)據(jù)向量排列組成一個(gè)新的mn×k階矩陣V，則V可以表示為

在非負(fù)矩陣分解的算法中，V＝WH＋ε。V為觀測圖像，W 為標(biāo)準(zhǔn)圖像，H 為權(quán)值（傳感器增益），ε為噪聲。圖像原始數(shù)據(jù)V可以近似分解為矩陣Wmn×r和矩陣Hr×k的乘積。故非負(fù)矩陣W 其實(shí)是對圖像原始數(shù)據(jù)V進(jìn)行線性逼近的一組基，H則是圖像原始數(shù)據(jù)V在W 上的非負(fù)投影系數(shù)。當(dāng)W的列數(shù)及特征基的數(shù)量r與實(shí)際數(shù)據(jù)集特征空間的維數(shù)一致時(shí)，得到的基向量是對圖像原始數(shù)據(jù)的最佳表示。這里考慮r＝1的特殊情況，此時(shí)唯一的特征基包含了圖像原始數(shù)據(jù)的整體特征，可用于源圖像的近似再現(xiàn)，而且可以有效抑制背景噪聲。

低頻子帶融合規(guī)則由dxJ（X＝A，B）構(gòu)造（6）式中矩陣V，然后對其進(jìn)行非負(fù)矩陣分解，分解后所得W 即為融合圖像的aFJ。

2.2 高頻子帶融合策略

圖像分解后高頻子帶中主要包含豐富的邊緣細(xì)節(jié)，如圖像的輪廓信息，目前尚未有精確的測度。若采用一般的直接選取絕對值最大或者區(qū)域特征值取大，都不同程度上造成圖像的不連續(xù)或者使得圖像的輪廓與非輪廓區(qū)域之間模糊，造成清晰度下降。這里采用高斯加權(quán)區(qū)域能量與區(qū)域標(biāo)準(zhǔn)差一致性選擇的融合策略，具體過程如下。

1）對于圖像 A、B 的高頻子帶系數(shù)｛HBj，k，HAj，k｝，以像素點(diǎn)（i，j）為中心，區(qū)域窗口為 M×N的加權(quán)區(qū)域能量及區(qū)域方差［7］的定義為

式中：X 取A 或B；μx（i，j）為區(qū)域均值；ω 表示區(qū)域加權(quán)權(quán)值矩陣，一般離像素中心點(diǎn)（i，j）越近，權(quán)值越大。

2）計(jì)算區(qū)域方差與加權(quán)區(qū)域能量相結(jié)合的融合指標(biāo)。對REXl，k和σXl，k先進(jìn)行歸一化處理，然后通過下式進(jìn)行計(jì)算：

式中：X取A 或B；a、b滿足a＋b＝1，一般取a＝b＝0.5。

3）確定融合后的系數(shù) HFj，k（i，j）：

3 實(shí)驗(yàn)分析

3.1 NSST變換的不同尺度和方向的分解

圖3 NSST變換的各層分解圖Fig.3 Each layer decomposition of NSST transform

非抽樣剪切波變換是剪切波變換（ST）的非下采樣變換方法，相比較ST而言，它具有平移不變性，穩(wěn)定性等特點(diǎn)。其對應(yīng)系數(shù)同樣能夠更準(zhǔn)確地描述圖像方向性和表面特征，更有效地提取源圖像中的有用信息［13］。下面以測試圖像（“zoneplate.png”）為例進(jìn)行2層尺度分解，方向參數(shù)L＝｛1，2｝，圖3為測試圖像的分解結(jié)果，對應(yīng)各層方向子帶個(gè)數(shù)分別為4、6。出的基于NSCT域融合算法——低頻采用區(qū)域相關(guān)系數(shù)和區(qū)域能量相結(jié)合，高頻采用加權(quán)區(qū)域方差的融合規(guī)則（NSCT－EN－WVAR）；ZHENG［7］等人提出的基于NSST域融合算法——低頻采用改進(jìn)平均梯度計(jì)算加權(quán)參量，高頻采用區(qū)域能量和區(qū)域方差結(jié)合的融合規(guī)則（NSST－GR－VAR）；本文融合算法此處記為NSST－OR。這里用視覺觀察對融合后的圖像進(jìn)行主觀評價(jià)，應(yīng)用信息熵EN、互信息MI和加權(quán)邊緣信息量QAB／F這3項(xiàng)指標(biāo)進(jìn)行客觀評價(jià)［12］。

下面將上述算法依次對紅外可見光圖像——圖3（a）和圖3（b）、多聚焦圖像——圖3（c）和圖3（d）進(jìn)行融合實(shí)驗(yàn)并進(jìn)行客觀評價(jià)。融合后的圖像依次為圖4（a）～圖4（f）以及圖5（a）～圖5（f）。首先從圖4的融合結(jié)果看出，DWT的融合效果在清晰度，邊緣有明顯的波紋干擾；基于曲波和輪廓波變換的效果圖像層次感略高，但仍然存在虛影模糊的情形；從圖4（d）和圖4（e）兩幅圖來看，圖中下方的樹木融合效果差異性較大，這和融合策略自身有一定的關(guān)系；圖4（f）為本文算法得到的融合結(jié)果，人物效果清晰，圖像整體對比度適中，細(xì)節(jié)信息也相對豐富。從圖5的融合結(jié)果可以看到：DWT的圖5（a）像整體比較模糊，而且邊緣細(xì)節(jié)信息不能較好地體現(xiàn)出來；圖5（b）和圖5（c）圖像清晰度上有一定的改善，但鐘表輪廓有明顯的偽吉布斯現(xiàn)象；圖5（d）～圖5（f）從視覺角度觀察圖像整體都比較清晰，邊緣刻畫都較為準(zhǔn)確，這說明在光譜差異不大的圖像融合過程中后面幾種算法都能得到較理想的結(jié)果。表1給出不同融合策略的客觀評價(jià)數(shù)據(jù)。

圖3 實(shí)驗(yàn)測試圖像Fig.3 Images of experiment test

3.2 結(jié)果分析

為驗(yàn)證本文所述融合算法的有效性，主要通過對嚴(yán)格配準(zhǔn)的紅外與可見光圖像以及多聚焦圖像進(jìn)行驗(yàn)證，首先選取傳統(tǒng)融合算法的情形下［10］，分別采取不同小波對進(jìn)行融合，即選取傳統(tǒng)融合策略（低頻子帶系數(shù)取平均，高頻系數(shù)取絕對值最大），分別在離散小波變換［9］（DWT－mean－max），曲 波 變 換［11］（curvelet－mean－max）、輪 廓 波 變 換（contourlet－mean－max）下進(jìn)行融合；其次針對不同融合策略進(jìn)行比較，這里選用ZHANG［8］等人提

圖4 紅外可見光不同融合策略的融合結(jié)果Fig.4 Infrared optical fusion results of different fusion strategies

圖5 多聚焦圖像不同融合策略的融合結(jié)果Fig.5 Multi－focus image fusion results of different fusion strategies

表1 融合結(jié)果的客觀評價(jià)數(shù)據(jù)Table 1 Objective evaluation data of fusion results

從表1結(jié)果可以看出，多聚焦圖像中曲波和輪廓波的傳統(tǒng)融合策略結(jié)果熵值雖然較高，但是圖像輪廓及邊緣信息的融合效果較差，而且主觀視覺來看明顯存在塊狀吉布斯效應(yīng)；而算法4和算法5的融合后圖像雖然邊緣信息關(guān)聯(lián)度有了明顯的提升，但包含平均信息量較低，綜合結(jié)果來看本文算法的融合效果較好。紅外與可見光圖像的各項(xiàng)指標(biāo)中，采用非下采樣的輪廓波，剪切波變換效果均優(yōu)于前3種融合效果，而本文算法在信息熵及圖像邊緣信息保留度上均有較為滿意的結(jié)果，進(jìn)一步驗(yàn)證了本文算法的有效性。

［1］ Do M N，Vetter M.The contourlet transform：an efficient directional multiresolution image repreisentation［J］.IEEE Trans.on Image Processing，2005，14（12）：2091－2106.

［2］ Cunha A L，Zhou J P，Dom N.The non－subsampled contourlet transform：theory，design，and application［J］.IEEE Trans.on Image Processing，2006，15（10）：3089－3101.

［3］ Guo K，Labate D.Characterization and analysis of edges using the continuous shearlet transform［J］.SIAM Journal on Imaging Sciences，2009，2（3）：959－986.

［4］ Cao Yuan，Li Shutao，Hu Jianwen.Multi－focus image fusion by nonsub sampled shear－let transform［C］.US：IEEE，2011：17－21.

［5］ Yang Yuetao，Zhu Ming，He Baigen，et.al.，F(xiàn)usion algorithm based on improved projected gradient NMF and NSCT［J］.Optics and Precision Engineering.2011，19（5）：1143－1145.

楊粵濤，朱明，賀柏根，等.采用改進(jìn)投影梯度非負(fù)矩陣分解和非采樣Contourlet變換的圖像融合方法［J］.光學(xué)精密工程，2011，19（5）：1143－1150.

［6］ Li Le，Zhang Yujin.A survey algorithms of nonnegative matrix factorization［J］.Acta Electronica Sinica，2008，36（4）：737－743.

李樂，章毓晉.非負(fù)矩陣算法綜述［J］.電子學(xué)報(bào)，2008，36（4）：737－743.

［7］ Zheng Wei，Sun Xueqing，Li Zhe，et al.Image fusion based on shearlet transform and region characteristics［J］.Laser Technology，2015，39（1）：50－56.

鄭偉，孫雪青，李哲，等.Shearlet變換和區(qū)域特性相結(jié)合的圖像融合［J］.激光技術(shù)，2015，39（1）：50－56.

［8］ Zhang Qiang，Guo Baolong.Fusion of infrared and visible light images based on nonsubsampled contourlet transform［J］.Journal of Infrared and Millimeter Waves，2007，26（6）：476－480.

張強(qiáng)，郭寶龍.一種基于非采樣Contourlet變換紅外圖像與可見光圖像融合算法［J］.紅外與毫米波學(xué)報(bào)，2007，26（6）：476－480.

［9］ Yue Zhen，Li Fanming.Polarization image fusion algorithm based on wavelet transform［J］.Journal of Applied Optics，2014，35（2）：321－326.

岳振，李范鳴.一種基于小波變換的紅外偏振融合［J］.應(yīng)用光學(xué)，2014，35（2）：321－326.

［10］Na Yan，Jiao Licheng.Image fusion method based multiresolution analysis theory［M］.Xi’an Press of Xidian University，2007.

那彥，焦李成.基于多分辨分析理論的圖像融合方法［M］.西安：西安電子科技大學(xué)出版社，2007.

［11］Alejaily A M，Rube I A E，Mangoud M A.Fusion of remote sensing images using contour let transform［M］.Berlin：Springer－Verlag，2008：213－218.

［12］Lu Yaning，Guo Lei，Li Huihui.Remote sensing image fusion using edge information and features of SAR image based on curve let transform［J］.Acta Photonic Sinica，2012，41（9）：1118－1123.

［13］Gao Shaoshu，Jin Weiqi，Wang Lingxue，et al.Objective quality assement of image fusion［J］.Journal of Applied Optics，2011，32（4）：671－677.

高紹姝，金偉其，王嶺雪，等.圖像融合質(zhì)量客觀評價(jià)方法 ［J］.應(yīng)用光學(xué)，2011，32（4）：671－677.