溫 佳，馬彩文，趙軍鎖，王彩玲

(1.中國科學(xué)院軟件研究所天基綜合信息系統(tǒng)重點(diǎn)實(shí)驗(yàn)室，100190北京;

2.中國科學(xué)院西安光學(xué)精密機(jī)械研究所，710119西安;3.西安石油大學(xué)計(jì)算機(jī)學(xué)院，710065西安)

干涉高光譜技術(shù)在航空航天遙感領(lǐng)域中是一種很有實(shí)用價(jià)值的技術(shù)，通過這種技術(shù)可以獲得觀測目標(biāo)的光譜信息與空間信息，因此在氣象、軍事、環(huán)境監(jiān)測和地質(zhì)等領(lǐng)域都有較廣泛的應(yīng)用.隨著干涉高光譜儀在中國“嫦娥”探月衛(wèi)星和近期環(huán)境探測衛(wèi)星中的成功搭載，干涉高光譜的數(shù)據(jù)壓縮成為近年來的研究焦點(diǎn)之一.

干涉高光譜圖像數(shù)據(jù)是一種三維圖像數(shù)據(jù)，其海量數(shù)據(jù)造成在數(shù)據(jù)存儲(chǔ)和有限帶寬信道上的傳輸有一定程度的困難，所以針對其數(shù)據(jù)本身特點(diǎn)設(shè)計(jì)出適用于干涉高光譜數(shù)據(jù)的高效數(shù)據(jù)壓縮方法十分必要.

近年來，越來越多的科學(xué)家針對干涉高光譜數(shù)據(jù)特性提出了有效的壓縮算法，基于譜間DPCM和整數(shù)小波的超光譜圖像無損壓縮在2008年提出［1］，針對干涉高光譜圖像幀間不穩(wěn)定移位的改進(jìn)變換算法在2011年提出［2］，針對干涉高光譜圖像幀間相關(guān)的自適應(yīng)光程差算法在2011 年提出［3］.

本文在文獻(xiàn)［2］基礎(chǔ)上，考慮到每幀 LASIS干涉高光譜數(shù)據(jù)均存在幅值較大的豎直干涉條紋，并且干涉條紋在每幀之間不存在移位，而除干涉條紋之外的背景圖像均存在方向性移位的特點(diǎn)，嘗試在幀序列方向采用不同方向自適應(yīng)得到最佳預(yù)測值，并且改變傳統(tǒng)三維提升小波的變換順序以消除干涉條紋所造成的冗余以達(dá)到更好的壓縮效果和壓縮性能.

1 干涉高光譜圖像的數(shù)據(jù)特性

干涉高光譜圖像數(shù)據(jù)具有與其他圖像不同的特性.圖1為LASIS干涉高光譜圖像三維示意圖.由于其特殊的推掃式成像原理，干涉高光譜圖像具有明顯的豎直干涉條紋存在.

圖1 LASⅠS干涉高光譜圖像數(shù)據(jù)

由圖1可以明顯看出干涉高光譜與一般的普通圖像存在明顯差異，干涉高光譜圖像具有如下特點(diǎn):

1)光的干涉是成像原理，干涉高光譜圖像并非光的直接成像，而是干涉圖像.圖像內(nèi)存在明顯的豎直干涉條紋，這些豎直的干涉條紋隨著光線調(diào)制程度的改變在不同位置呈現(xiàn)不同的強(qiáng)度.

2)干涉高光譜圖像是三維數(shù)據(jù)，具有多維方向的相關(guān)性.在聯(lián)合調(diào)制干涉高光譜圖像中，目標(biāo)是一個(gè)整體的推掃平移過程，幀內(nèi)圖像具有空間信息，在圖像幀間存在著平移現(xiàn)象(干涉條紋不存在).

2 提升小波變換

Haar變換的提升版本是最簡單的提升小波變換之一，稱為 S 變換［4］.

S變換步驟如下:

式中:x［n］為輸入信號(hào)，s［n］為提升小波變換后的低頻信號(hào)，d［n］為提升小波變換后的高頻信號(hào)，?」表示向下取整運(yùn)算．經(jīng)過S提升后的低頻信號(hào)與原始輸入信號(hào)的動(dòng)態(tài)范圍相同，高頻信號(hào)將變得很小，Said等［5］之后又提出了S+P(S變換 +Predication)，Zandi等［6］提出的 TS變換是基于梯形結(jié)構(gòu)的比S變換具有更高消失矩的整數(shù)小波變換.Sweldends等［7-11］在 Donoho 等［12］的基礎(chǔ)上又提出了采用提升框架的第二代小波變換結(jié)構(gòu)之后，Calderbank等［13］提出了基于提升格式的整數(shù)小波變換.

提升小波變換分為3個(gè)步驟:

1)lazy變換;

2)預(yù)測 (predication);

3)更新 (update).

第一步lazy變換是把數(shù)據(jù)分為偶數(shù)集合和奇數(shù)集合

第二步預(yù)測是用偶數(shù)集來預(yù)測奇數(shù)集，將產(chǎn)生的誤差作為高頻系數(shù).

第三步更新是用這些高頻系數(shù)來更新偶數(shù)集合作為低通系數(shù).

式中:si［n］為第i層原始提升的低通分量，di［n］為第i層對偶提升的高通分量．

3 根據(jù)干涉高光譜特性改進(jìn)的自適應(yīng)方向角三維提升小波變換

由于干涉多光譜圖像幀與幀之間有很強(qiáng)的方向性，并且三維LASIS圖像序列幀與幀之間的平移量是不穩(wěn)定的，實(shí)驗(yàn)數(shù)據(jù)表明一般相鄰幀之間的平移在1～3列不等.LASIS通常搭載在飛行器上進(jìn)行探測，運(yùn)用大孔徑靜態(tài)干涉成像光譜儀面陣探測器依靠推掃獲得二維的空間信息與一維的光譜信息，形成三維立體圖像數(shù)據(jù).在干涉高光譜圖像中具有明顯的豎直干涉條紋，且豎直干涉條紋的位置固定，不存在移位現(xiàn)象，豎直干涉條紋中含有光譜信息，并且干涉高光譜圖像幀間具有明顯平移，方向與探測器推掃方向基本一致，然而由于非勻速的推掃，平移方向存在抖動(dòng).

本文采用傳統(tǒng)5/3提升小波，在幀內(nèi)的行方向和列方向直接采用傳統(tǒng)5/3提升小波進(jìn)行小波變換:

為了提高運(yùn)算速度，對于r×c×f的LASIS三維數(shù)據(jù)，在列方向的提升小波變換之前，先將圖像拼接成［r，c×f］的2維矩陣，拼接方式如圖2所示.新生成的2維矩陣有c×f列，每列含有r個(gè)像素值，將新生成的［r，c×f］二維矩陣對應(yīng)式(2)中的x變量，將行數(shù)r對應(yīng)公式(2)的變量n，之后按照式(2)生成s0［n］，d0［n］的二維矩陣如圖 3 所示．將生成的 s0［n］、d0［n］按照式(5)、(6)進(jìn)行一次一維垂直方向的提升小波變換，結(jié)果如圖4所示，之后再把新生成的矩陣還原成原來的［r，c，f］的三維矩陣，即可完成對 r× c×f的LASIS三維數(shù)據(jù)的列方向提升小波變換.水平方向提升小波變換按照與之類似流程完成.對于幀序列方向的提升小波變換，則把三維LASIS數(shù)據(jù)拼接成［f，r×c］的2維矩陣，每k行的r×c個(gè)像素值是原始三維矩陣的第k幀的幀內(nèi)r×c個(gè)像素按水平方向拼接而成，k=1，2…f，由于LASIS數(shù)據(jù)幀間的方向特性，拼接后的2維［f，r×c］矩陣如圖5所示.新生成的2維矩陣有r×c列，每列含有f個(gè)像素值，按式(2)將新生成的［f，r×c］二維矩陣對應(yīng)式(2)中的x變量，將行數(shù)f對應(yīng)式(2)的變量n，之后按照公式(2)生成的s0［n］、d0［n］為的二維矩陣如圖6所示．將生成的 s0［n］、d0［n］按照式(5)、(6)進(jìn)行一次一維垂直方向的提升小波變換，結(jié)果如圖7所示．

圖2 垂直方向提升小波變換LSSIS 數(shù)據(jù)拼接示意

圖3 垂直方向提升小波變換中s0［n］和d0［n］示意

圖4 垂直提升小波變換后LASⅠS數(shù)據(jù)示意

圖5 幀間提升小波變換LASⅠS數(shù)據(jù)拼接示意

圖6 幀間方向提升小波變換中s0［n］和d0［n］示意

圖7 幀間提升小波變換后LASⅠS數(shù)據(jù)示意

之后再把新生成的矩陣還原成原來的［r，c，f］的三維矩陣，即可完成對r× c× f的 LASIS三維數(shù)據(jù)的幀序列方向提升小波變換.

由于干涉高光譜圖像的成像原理所導(dǎo)致其自身的特殊性質(zhì)，對于幀序列方向，由于干涉高光譜圖像數(shù)據(jù)在幀間存在1～3的移位，這在圖5中可以明顯看到，況且在進(jìn)行完一級(jí)提升小波變換之后，再做下一級(jí)提升小波變換時(shí)，幀間的移位會(huì)有加倍的變化，這一點(diǎn)在圖5、7的對比中也可以明顯看到.

直接對上述矩陣在垂直方向做提升小波變換在第二步的預(yù)測步驟并不能得到很好的結(jié)果，所得的高頻系數(shù)會(huì)偏大.本文針對LASIS三維數(shù)據(jù)的幀序列方向?qū)ｉT提出一種自適應(yīng)預(yù)測方向角的提升小波變換以適應(yīng)LASIS幀間特有的方向性特點(diǎn).

在預(yù)測步驟采用5/3的提升小波變換d1［n］=，d0［n］實(shí)際上是夾在 s0［n］、s0［n+1］中間的層，由于幀間1 ～ 3移位的存在，體現(xiàn)在d0［n］、s0［n］、s0［n+1］上為:s0［n］為了與d0［n］在垂直方向?qū)?yīng)，應(yīng)向左移動(dòng)1 ～3位，s0［n+1］為了與d0［n］在垂直方向?qū)?yīng)，應(yīng)向右移動(dòng)1～3位．

s0［n］和 s0［n+1］與 d0［n］最佳匹配的移位各有4種可能，其兩兩組合共有16種可能性存在，這樣進(jìn)行預(yù)測步驟會(huì)得到16種可能的d1i［n］，i=1，2，…，16，將d1i［n］的每一行的r× c個(gè)像素值的絕對值相加，會(huì)得到16個(gè)［f/2，1］的列，將其轉(zhuǎn)置成行并用一個(gè)［16，f/2］的二維矩陣d1_sum存儲(chǔ)，d1_sum有f/2列，每一列含有16種可能的值，保存每列最小值的索引值，f/2列會(huì)產(chǎn)生f/2個(gè)index索引值，index(k)∈［1，16］，k=1，2，…，f/2．如此可以快速的得到每幀的上一幀和下一幀對應(yīng)的移位最佳組合，新生成的d1［n］的每k行為 d1index(k)的第k行，k=1，2，…，f/2．在信道傳輸數(shù)據(jù)時(shí)，索引值需要與編碼后的碼字一起傳輸．

LASIS干涉高光譜三維數(shù)據(jù)由于豎直干涉條紋的存在并且有較大的幅值，如果為了提高運(yùn)算速度把提升小波的預(yù)測步驟中的 d0［n］、s0［n］作為整體進(jìn)行運(yùn)算，會(huì)影響到幀間的方向性，因?yàn)楸M管大部分背景像素存在整體的方向性移位，但是豎直的干涉條紋不存在移位并且具有較大的幅值.因此，在這里不采用傳統(tǒng)的對稱三維小波變換，先對LASIS數(shù)據(jù)進(jìn)行3級(jí)的一維垂直方向的提升小波變換，可以基本上消除掉LASIS數(shù)據(jù)中特有的豎直干涉條紋所帶來的影響，圖8所示為對16幀LASIS數(shù)據(jù)進(jìn)行3級(jí)豎直方向提升小波變換后的3幀小波系數(shù)示意圖.

由圖8可見，在低頻分量依然存在不移位的干涉條紋，而高頻分量的干涉條紋冗余已基本消除，由于本文采用拼接后整體進(jìn)行運(yùn)算，如果直接對3級(jí)垂直方向提升小波變換后的系數(shù)進(jìn)行幀序列方向預(yù)測，并不能得到最佳的移位索引.因此本文對含有干涉條紋冗余的低頻分量和干涉條紋冗余以基本被消除的高頻分量分別進(jìn)行3級(jí)自適應(yīng)幀序列方向的提升小波變換，這樣會(huì)分別得到存在干涉條紋影響的最佳移位索引和不存在干涉條紋影響的最佳移位索引，最后再將所得的結(jié)果進(jìn)行3級(jí)水平方向的提升小波變換.

本文提出自適應(yīng)預(yù)測提升小波變換算法的流程圖如圖9所示.

圖8 三級(jí)垂直方向提升小波變換后小波系數(shù)示意圖

圖9 本文提出算法的流程圖

4 實(shí)驗(yàn)結(jié)果

對兩組512×288×16的12 bvit LASIS三維數(shù)據(jù)在bpppb(比特/像素/波段)=0.1～0.3進(jìn)行實(shí)驗(yàn)，結(jié)果如表1～2所示.

表1 不同bpppb下LASⅠS數(shù)據(jù)1的信噪比 dB

表2 不同bpppb下LASⅠS數(shù)據(jù)2的信噪比 dB

由表1、2可知，使用本文提出的變換方法相對于傳統(tǒng)提升小波變換方法，在經(jīng)過3DSPECK編碼后，在給定編碼比特?cái)?shù)目即壓縮比一定的情況下重構(gòu)圖像可以獲得更高的信噪比.圖10、11顯示了1組512×288×16的LASIS干涉高光譜圖像數(shù)據(jù)經(jīng)過本文提出改進(jìn)的變換方法與傳統(tǒng)三維對稱提升小波變換和2011年提出的自適應(yīng)小波變換［2］后經(jīng)過3DSPECK編碼后重構(gòu)結(jié)果的比較.將3種方法得到的LASIS重構(gòu)干涉高光譜數(shù)據(jù)幀通過后期處理得到LSMIS(Large Spatially Modulated Interference Spectral Image)干涉圖后，還原光譜曲線如圖12所示，還原光譜曲線與原始光譜曲線的均方誤差如表3所示.

圖10 LASⅠS原始數(shù)據(jù)

圖11 比特率為0.3時(shí)3種不同方法的恢復(fù)圖像

圖12 3種不同方法的恢復(fù)光譜曲線

表3 不同方法下還原光譜曲線的均方誤差

5 結(jié)語

從本文實(shí)驗(yàn)結(jié)果可以明顯看出如果仍然運(yùn)用傳統(tǒng)小波變換順序進(jìn)行幀序列方向的自適應(yīng)方向角預(yù)測，對高頻小波系數(shù)間方差減小并不明顯，因?yàn)楦缮鏃l紋的存在影響了整體自適應(yīng)方向預(yù)測的結(jié)果，通過本文方法消除垂直方向冗余并且對含有干涉條紋冗余的低頻子帶和不含干涉條紋冗余的高頻子帶分別進(jìn)行自適應(yīng)方向預(yù)測后，高頻子帶的小波系數(shù)間方差得到了大幅度改善.將提出的小波變換方法與傳統(tǒng)方法產(chǎn)生的小波系數(shù)在指定碼率下進(jìn)行3DSPECK編碼，大量實(shí)驗(yàn)數(shù)據(jù)表明本文提出的改進(jìn)方法的重構(gòu)圖像可以獲得更高的信噪比，且后期處理得到的光譜曲線相對于原始光譜曲線具有更小的均方誤差.

［1］吳冬梅，王軍，張海寧.基于譜間 DPCM和整數(shù)小波的超光譜圖像無損壓縮［J］.光子學(xué)報(bào)，2008，37(1):156-159.

［2］WEN J，MA C W，SHUI P L.A 3D non-linear orientation prediction wavelet transform for interference hyperspectral images compression［J］.Optics Communications，2011，284(7):1770-1777.

［3］WEN J，MA C W，SHUI P L.An adaptive OPD and dislocation prediction used characteristic of interference pattern for interference hyperspectral image compression［J］.Optics Communications，2011，284(20):4903-4909.

［4］BLUME H，F(xiàn)AND A.Reversible and irreversible image data compression using the S.transform and Lempel-Ziv coding［C］//Proceedings of SPIE Conference on Medical Imaging.Newport Beach:CA，1989，1091:2.

［5］SAID A，PEARLMAN W A.An image multiresolution representation for lossless and lossy compression［J］.IEEE Transactions on Image Processing，1996，5(9):1303-1310.

［6］ZANDI A，ALLEN J，SCHWARTZ J E.Compression with reversible embedded wavelets［C］//Proeeedings of IEEE Data Compression Conference.Snowbird:DCC，1995:212

［7］SWELDENS W.The lifting seheme:a new philosophy in biorthogonal wavelet constructions［C］//Proceedings of SPIE Conference on Wavelet Applications in Signal and Image Processing.San Diego:CA，1995，2569:68.

［8］SWELDENS W，SHROEDER P.Building your own wavelets at home［R］.Columbia:University of South Carolina，1995.

［9］SWELDENS W.The lifting scheme:a custom-design construction of biothogonal wavelets［J］.Applied and Computational Harmonic Analysis，1996，3(2):186-200.

［10］SWELDENS W.The lifting scheme:construction of second generation wavelets［J］. SLAM Jounalof Mathematical Analysis，1997，29(2):511-546.

［11］DAUBEEHIES I，SWELDENS W.Factoring wavelet transforms into lifting steps［J］.Joumal of Fourier Analysis and Applications，1998，4(3):247-269.

［12］DONOHO D L.Interpolating wavelet transform［M］.PaloAltoCounty:DepartmentofStatistics，Standford University，1992.

［13］CALDERBANK A R，DAUBECHIES I，SWELDENS W.Wavelet transforms that map integers to integers［J］.Applied and Computational Harmonic Analysis，1998，31(5):332-369.