吳曉萍,楊武年*,李國明
(1.成都理工大學(xué) 地學(xué)空間信息技術(shù)國土資源部重點實驗室/遙感與GIS研究所,成都 610059;2.四川省第三測繪工程院,成都 610500)
隨著空間技術(shù)的發(fā)展,從不同物理特性的傳感器獲得的海量不同空間分辨率、光譜分辨率和時間分辨率的遙感影像,形成了多級多分辨率的影像金字塔序列,實時的為用戶提供對地觀測數(shù)據(jù)源。每種傳感器所獲得的遙感數(shù)據(jù)只能從一個或幾個方面來反映事物的特性,而遙感影像的信息提取及應(yīng)用通常要求把多傳感器、多光譜和多分辨率影像結(jié)臺起來分析,以克服遙感影像在解譯過程中單一信息源不足的問題[1]。采取有效方法,對多光譜遙感數(shù)據(jù)進(jìn)行融合,使圖像同時具有較高的光譜和空間分辨率,提高了圖像的視覺效果和圖像特征識別及分類的精度。
ZY-3是2012年1月9日中國發(fā)射的首顆民用高分辨率光學(xué)傳輸型立體測圖衛(wèi)星[2]。主要搭載有一臺地面分辨率2.1m的高分辨率正視全色延時積分成像(TDI CCD)相機(jī)、兩臺地面分辨率3.6m的前視、后視全色(TDI CCD)相機(jī)和一臺地面分辨率5.8m的正視多光譜相機(jī),其中多光譜通道分別獲取波長為0.45μm~0.52μm藍(lán)光光譜信息、0.52μm~0.59μm綠光光譜信息、0.63μm~0.69 μm紅光光譜信息和0.77μm~0.89μm近紅外光譜信息。
目前,國內(nèi)、外許多學(xué)者應(yīng)用HIS變換、乘積變換、主成份變換、小波變換等不同的融合方法,對遙感數(shù)據(jù)融合技術(shù)進(jìn)行了研究,對于SPOT、IKO-NOS、Quick Bird、Landsat、ETM+等影像參與融合的研究較多,尚無對ZY-3數(shù)據(jù)融合的研究。因此,作者采用Brovey變換、主成分變換、IHS變換、小波變換、Gram-Schimdt光譜銳化等五種目前使用較為廣泛的圖像融合技術(shù),將資源三號的空間分辨率5.8m多光譜數(shù)據(jù)與空間分辨率2.1m的正視全色數(shù)據(jù)進(jìn)行融合,提高影像分辨率。通過研究可行的資源三號數(shù)據(jù)融合方法,使其發(fā)揮更大的實用價值。
Brovey圖像融合是一種通過歸一化后的多光譜影像與高分辨率影像乘積來增強(qiáng)影像空間信息的融合方法。其定義式如公式(1)所示[3]:
式中 DNF為融合后影像像素值;DNLi(i=1,2,3)為多光譜影像中的第i波段的像素值;DNH為全色影像波段的像素值。
主成分變換是將多個變量通過線性變換選出幾個重要變量來進(jìn)行多元統(tǒng)計分析的一種方法,并在其特征基礎(chǔ)上進(jìn)行的多波段影像的正交線性變換[4]。正反變換式如公式(2)、公式(3)所示。
式中 LP為多光譜影像L經(jīng)過主成分正變換得到的主成分影像,T為主成分正變換矩陣。
式中 H為融合后的影像,T-1為反變換矩陣,是正變換矩陣的逆矩陣。
在色度學(xué)中,把彩色影像的紅(R)、綠(G)、藍(lán)(B)變換成明度(I)、色度(H)、飽和度(S)稱為正變換,反之為反變換。其正變換式如公式(4)-公式(6)所示[4]。
其反變換式如公式(7)所示[4]:
式中 I1、I2為中間變量。根據(jù)公式(4)-公式(6)可知,先將影像進(jìn)行RGB系統(tǒng)到IHS系統(tǒng)的空間變換;再將全色波段數(shù)據(jù)與分解出的I分量進(jìn)行直方圖匹配,用其匹配結(jié)果全色波段替代I分量與H分量、S分量合成;最后根據(jù)公式(7)將其反變換到RGB空間中 。
將多光譜遙感影像和全色光遙感影像在不同的頻域率進(jìn)行多次分解,然后采用融合算法對小波系數(shù)進(jìn)行融合,得到包含多光譜信息與全色光信息的小波系數(shù),最后對此小波系數(shù)進(jìn)行逆變換,得到融合圖像。最常用的融合方法是小波替代法,它是將多光譜圖像的細(xì)節(jié)用全色圖像的細(xì)節(jié)來代替,再將剩余的多光譜影像進(jìn)行融合,其計算式如公式(8)公式(10)所示[3]:
Gram-Schimdt光譜銳化是線性多元統(tǒng)計和代數(shù)中常用的方法,它是通過多維影像進(jìn)行正交變換來消除冗余的信息,利用多光譜影像模擬全色高分辨率影像,并將其作為Gram-Schimdt變換的第一個分量,通過調(diào)整高分辨率全色影像的統(tǒng)計值來匹配Gram-Schimdt變換后的第一個分量,再對其進(jìn)行逆變換,變換后的影像各分量的信息并沒有發(fā)生變化,其計算公式如公式(11)所示[3]。
式中 GST是弦GS變換后產(chǎn)生的第T個分量;ΒΤ是原始多光譜影像的第T個波段;uT是第T個原始多光譜波段的灰度均值;?(BT,GSl)是協(xié)方差。
GS逆變換的公式如公式(12)所示:
圖1 原始多光譜影像圖Fig.1 Original multispectral image
圖2 區(qū)域放大全色影像Fig.2 Partial enlarged panchromatic image
圖3 區(qū)域放大多光譜影像Fig.3 Partial enlarged multispectral image
利用Brovey融合方法、主成分變換、IHS變換、小波變換、Gram-Schimdt光譜銳化法等五種方法,以四川省某區(qū)域資源三號遙感影像為數(shù)據(jù)源進(jìn)行了融合實驗。為了比較不同融合方法的融合效果,從整個研究區(qū)域中選擇一個區(qū)域,并以1∶1比例尺顯示。圖1為示范區(qū)資源三號多光譜影像。圖2為區(qū)域放大的全色影像,分辨率為2.1m。圖3為區(qū)域放大的多光譜影像,分辨率為5.8m。圖4圖8分別是采用上述五種融合方法處理后的區(qū)域放大影像。
根據(jù)美國查維茨Chavez等(1984)提出的最佳指數(shù)(OIF)概念,選擇OIF最大的波段組合,可以得到3波段、2波段、1波段組合的信息量最大。因此評價所用影像均為3(R)、2(G)、1(B)波段真彩色合成。從圖3-圖8影像分辨率上可以看出:①經(jīng)過融合處理后的影像空間分辨率有了很大地提高,房屋、道路、農(nóng)田的界線等都清晰可見;②從提高空間紋理信息及影像亮度上來看,幾種融合方法并無明顯差別,比原始多光譜影像有了很大提高;③從色調(diào)上看,除brovey融合方法外,其他融合方法所得到的影像與原始影像色調(diào)相當(dāng),保留了其原有色調(diào)信息。
圖4 Brovey融合影像Fig.4 Brovey fusion image
圖5 主成分變換影像Fig.5 The principal component transform image
圖6 IHS變換影像Fig.6 IHS transform image
圖7 小波變換影像Fig.7 Wavelet transform image
根據(jù)圖像融合前、后目視判別對比作出定性評價是最簡單、最直接的評價方法,但不同融合方法產(chǎn)生的光譜失真可能會導(dǎo)致不可靠的判別和應(yīng)用。因此為進(jìn)一步客觀定量地評價融合效果,作者將綜合利用兩類統(tǒng)計參數(shù)對融合效果進(jìn)行評價:①反映影像空間細(xì)節(jié)特性參數(shù)(標(biāo)準(zhǔn)差、信息熵和平均梯度);②反映光譜保持特性參數(shù)(偏差指數(shù)、相關(guān)系數(shù)和光譜扭曲程度)。
2.3.1 標(biāo)準(zhǔn)差
標(biāo)準(zhǔn)差反映了圖像中每個像元的灰度值與灰度平均值的離散程度,標(biāo)準(zhǔn)差越大說明圖像中灰度值分布越分散,影像的灰度直方圖越趨于水平,圖像所含的信息量越豐富。在數(shù)理統(tǒng)計中,標(biāo)準(zhǔn)差可用公式(13)計算[5]:
圖8 GS光譜銳化影像Fig.8 GS spectral sharpening image
式中 k為標(biāo)準(zhǔn)差;M、N 為圖像的行列數(shù);h(i,j)為第i行、第j列對應(yīng)像元的灰度值;ˉh為圖像平均灰度值。計算過程由IDL編程實現(xiàn),其結(jié)果見表1。
2.3.2 信息熵
信息熵是用來衡量融合影像信息量豐富程度的重要指標(biāo),在影像的表示上為偏離影像直方圖高峰值灰度區(qū)的大小,通過對圖像信息熵的比較可以對比出圖像的細(xì)節(jié)表現(xiàn)能力,融合影像的熵值越大其信息量也就越大。對于灰度范圍{0,1,…,N-1}的圖像直方圖,其熵的定義式如公式(14)所示[6]:
表1 融合前后影像標(biāo)準(zhǔn)差對比Tab.1 Standard deviation contrast before and after the fusion image
式中 pi為第i個灰度的出現(xiàn)概率。計算過程由Matlab編程實現(xiàn),其結(jié)果見表2。
表2 融合前后影像信息熵對比Tab.2 Information entropy contrast before and after the fusion image
2.3.3 平均梯度
影像質(zhì)量的改進(jìn)用平均梯度來評價,梯度越大,則影像的清晰程度越高,邊緣細(xì)節(jié)也就更加清晰,其計算式如公式(15)所示[6]:
式中 M、N 為圖像的行列數(shù);h(i,j)為第i行、第j列對應(yīng)像元的灰度值。計算過程由IDL編程實現(xiàn),其結(jié)果見表3。
2.3.4 偏差指數(shù)
偏差指數(shù)反映了影像融合前、后對應(yīng)多光譜波段之間的光譜質(zhì)量差異,是融合后的影像與多光譜影像灰度值差值的絕對值與多光譜影像的灰度值之比,比值越小說明融合影像與原多光譜影像的偏離程度越小,光譜信息繼承越好,光譜特征就得到了最大限度的保留,其計算式如公式(16)所示[7]:式中 M、N 為遙感影像的行列數(shù);Ia(i,j)、Ib(i,j)分別為融合前、后同一波段相同位置第i行、第j列對應(yīng)像元的灰度值。計算過程由IDL編程實現(xiàn),其結(jié)果見表4。
表3 融合前后影像平均梯度對比Tab.3 Average gradient contrast before and after the fusion image
表4 融合前后影像偏差指數(shù)對比Tab.4 Deviation index contrast before and after the fusion image
2.3.5 相關(guān)系數(shù)
圖像的相關(guān)系數(shù)反映了兩幅圖像的相關(guān)程度,通過比較融合前、后的圖像相關(guān)系數(shù),就可以看出多光譜圖像的光譜信息的改變程度,相關(guān)系數(shù)越大,越能保持原多光譜影像的光譜特征。兩幅圖像的相關(guān)系數(shù)定義式如公式(17)所示[8-9]:
式中 M、N 分別為圖像行列數(shù);Ia(i,j)、Ib(i,j)分別為融合前、后同一波段相同位置第i行、第j列對應(yīng)像元的灰度值分別為原圖像和融合后圖像的均值。計算過程由IDL編程實現(xiàn),其結(jié)果見表5。
2.3.6 光譜扭曲程度
圖像光譜扭曲程度直接反映了多光譜圖像的光譜失真程度,光譜扭曲定義式如公式(18)所示[10]:
式中 R 為光譜扭曲值;n為圖像的大??;Ia(i,j)、Ib(i,j)分別為融合前、后同一波段相同位置第i行、第j列對應(yīng)像元的灰度值。計算過程由IDL編程實現(xiàn),其結(jié)果見表6。
表5 融合前后影像相關(guān)系數(shù)對比Tab.5 Correlation coefficient contrast before and after the fusion image
表6 融合前后影像光譜扭曲程度對比Tab.6 The contrast of spectral distortion degree before and after the fusion image
通過圖3-圖8、表1-表6中五種融合結(jié)果與原始多光譜影像對比可知,五種融合方法均能較大程度提高多光譜圖像的清晰度,提高數(shù)據(jù)的空間分辨率。但通過定量評價表明:①Brovey融合影像色調(diào)偏離較大,所含信息量少、光譜失真過大,圖像邊緣信息減弱(圖4、表3、表4、表6);②IHS融合影像光譜失真過大,圖像邊緣信息減弱(圖6、表3、表6);③小波變換融合影像空間分辨率和信息量兩方面改善效果較差(圖7、表1、表2);④GS光譜銳化和主成分變換融合的影像從幾個評價值中可以看出(除平均梯度外),都明顯優(yōu)于其他的數(shù)據(jù)融合方法。主成分分析算法將高分辨率圖像替換第1主成分,保留了大量的全色圖像信息,該算法雖然提高了多光譜圖像的分辨率和平均梯度(圖5、表3),但是由于舍棄了過多的原多光譜數(shù)據(jù)信息,導(dǎo)致光譜失真程度比GS光譜銳化融合影像偏大,影像相關(guān)程度偏低(表4-表6)。
作者采用多種融合方法,對資源三號正視全色和多光譜影像進(jìn)行了融合試驗。從主觀定性和客觀定量兩個方面來分析、評價了五種融合結(jié)果,由上綜合分析可看出:在資源三號正視全色與多光譜影像融合時,Gram-Schimdt光譜銳化方法綜合效果優(yōu)于其它方法,融合的圖像信息量更加豐富,且光譜特征保留較好,主成分變換方法次之。這對遙感影像解譯、分類精度的提高,專題圖件的制作、有效地利用不同尺度、不同時相、不同類型的海量遙感數(shù)據(jù)及拓展資源三號影像的進(jìn)一步應(yīng)用范圍,具有重要的現(xiàn)實意義。
[1]賈永紅,李德仁,孫家柄.多源遙感影像數(shù)據(jù)融合[J].遙感技術(shù)與應(yīng)用,2000,15(1):41-44.
[2]李德仁.我國第一顆民用三線陣立體測圖衛(wèi)星-資源三號測繪衛(wèi)星[J].測繪學(xué)報,2012,41(3):317-322.
[3]章孝燦,黃智才,趙元洪.遙感數(shù)字圖像處理[M].杭州:浙江大學(xué)出版社,1997.
[4]MARFA GONZALEZ-AUDICANA, JOSE LUIS SALETA,RAQUEL GACAIA.Fusion of multispectral andpanchromatic images using improved IHS and PCA merges based on wavelet decomposition[J].IEEE Transactions on Geoscience and Remote Sensing,2004,42(6):1291-1298.
[5]劉貴喜,楊萬海.基于多尺度對比塔的圖像融合方法及性能評價[J].光學(xué)學(xué)報,2001,21(11):112-114.
[6]張寧玉,吳泉源.Brovey融合與小波融合對QuickBird圖像的信息量影響[J].遙感技術(shù)與應(yīng)用,2006,21(1):67-70.
[7]陳德超,周海波,陳中原,等.TM與SPOT影像融合算法比較研究[J].遙感技術(shù)與應(yīng)用,2001,16(2):110-115.
[8]李俊杰,李杏朝,傅俏燕.CBERS-02B星 HR與多光譜影像融合及評價[J].國土資源遙感,2007(2):43-47.
[9]LI G Q,LI X P,WANG H,et al.Cross-comparison between China HJ1A-CCD and Landsat TM data[J].2011IEEE International Geoscience and Remote Sensing Symposium Vancouver,Canada.23-29July,2011:4130-4133.
[10]許榕峰,徐涵秋.ETM+全色波段及其多光譜波段圖像的融合應(yīng)用[J].地球信息科學(xué),2004,6(1):99-103.