自適應提升小波在干涉高光譜壓縮中的應用

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

(1.中國科學院軟件研究所天基綜合信息系統重點實驗室，100190北京;

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

干涉高光譜技術在航空航天遙感領域中是一種很有實用價值的技術，通過這種技術可以獲得觀測目標的光譜信息與空間信息，因此在氣象、軍事、環境監測和地質等領域都有較廣泛的應用.隨著干涉高光譜儀在中國“嫦娥”探月衛星和近期環境探測衛星中的成功搭載，干涉高光譜的數據壓縮成為近年來的研究焦點之一.

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

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

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

1 干涉高光譜圖像的數據特性

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

圖1 LASⅠS干涉高光譜圖像數據

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

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

2)干涉高光譜圖像是三維數據，具有多維方向的相關性.在聯合調制干涉高光譜圖像中，目標是一個整體的推掃平移過程，幀內圖像具有空間信息，在圖像幀間存在著平移現象(干涉條紋不存在).

2 提升小波變換

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

S變換步驟如下:

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

提升小波變換分為3個步驟:

1)lazy變換;

2)預測 (predication);

3)更新 (update).

第一步lazy變換是把數據分為偶數集合和奇數集合

第二步預測是用偶數集來預測奇數集，將產生的誤差作為高頻系數.

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

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

3 根據干涉高光譜特性改進的自適應方向角三維提升小波變換

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

本文采用傳統5/3提升小波，在幀內的行方向和列方向直接采用傳統5/3提升小波進行小波變換:

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

圖2 垂直方向提升小波變換LSSIS 數據拼接示意

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

圖4 垂直提升小波變換后LASⅠS數據示意

圖5 幀間提升小波變換LASⅠS數據拼接示意

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

圖7 幀間提升小波變換后LASⅠS數據示意

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

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

直接對上述矩陣在垂直方向做提升小波變換在第二步的預測步驟并不能得到很好的結果，所得的高頻系數會偏大.本文針對LASIS三維數據的幀序列方向專門提出一種自適應預測方向角的提升小波變換以適應LASIS幀間特有的方向性特點.

在預測步驟采用5/3的提升小波變換d1［n］=，d0［n］實際上是夾在 s0［n］、s0［n+1］中間的層，由于幀間1 ～ 3移位的存在，體現在d0［n］、s0［n］、s0［n+1］上為:s0［n］為了與d0［n］在垂直方向對應，應向左移動1 ～3位，s0［n+1］為了與d0［n］在垂直方向對應，應向右移動1～3位．

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

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

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

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

圖8 三級垂直方向提升小波變換后小波系數示意圖

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

4 實驗結果

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

表1 不同bpppb下LASⅠS數據1的信噪比 dB

表2 不同bpppb下LASⅠS數據2的信噪比 dB

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

圖10 LASⅠS原始數據

圖11 比特率為0.3時3種不同方法的恢復圖像

圖12 3種不同方法的恢復光譜曲線

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

5 結語

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

［1］吳冬梅，王軍，張海寧.基于譜間 DPCM和整數小波的超光譜圖像無損壓縮［J］.光子學報，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，FAND 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.