基于全變分正則項的CASSI數據重構算法

王業超 陳曉麗 鐘曉明 趙海博 張麗莎 蘇云

基于全變分正則項的CASSI數據重構算法

王業超1,2陳曉麗1,2鐘曉明1,2趙海博1,2張麗莎1,2蘇云1,2

（1 北京空間機電研究所，北京 100094） （2 先進光學遙感技術北京市重點試驗室，北京 100094）

在編碼孔徑快照光譜儀（Coded Aperture Snapshot Spectral Imager，CASSI）光譜重構算法設計中，兩步迭代收縮閾值（Two-Step Iterative Shrinkage/Threshold，TwIST）算法實現了CASSI的光譜重構，但因正則項選取和引入的噪聲隨迭代次數而不斷放大，導致其光譜重構精度低于80%。文章在TwIST算法基礎上，以光譜圖像具有空間平滑過渡特性為先驗知識，提出兩點改進：一是選擇全變分正則約束項；二是對每一步迭代的更新項進行全變分去噪處理。為了驗證改進后的算法，文章通過計算機仿真CASSI的采集數據，得出仿真數據光譜重構精度為90.93%；并根據CASSI樣機采集試驗數據，得出試驗數據的光譜重構精度為86.56%。改進后的算法可為以后CASSI數據重構提供參考。

兩步迭代收縮閾值 壓縮感知 編碼孔徑 全變分 光譜重構 光譜儀

0 引言

編碼孔徑快照光譜儀（Coded Aperture Snapshot Spectral Imager，CASSI）是在壓縮感知理論框架下發展起來的。壓縮感知理論指出若信號可壓縮或在某個變換域稀疏，則利用與變換基不相關的觀測矩陣將高維信號投影至低維空間，可通過求解優化問題從少量投影信息中以高概率重構出原始信號，并證明投影信息包含重構所需的足夠信息[1-5]。在壓縮感知理論框架下，只要信號滿足可壓縮或稀疏的條件，即可通過遠低于奈奎斯特采樣頻率采樣，重構出原始信號[6-8]。上述兩個結論為CASSI算法設計提供理論支撐?；趬嚎s感知理論，文獻[9-13]提出并發展CASSI及其光譜重構算法的設計，但光譜重構精度較低。CASSI系統核心要素和難點是從探測器采樣的二維混疊圖中重構三維光譜數據立方體。文獻[14-20]提出通過硬件與算法聯合設計的思路重構CASSI光譜數據。受聯合設計啟發，研究人員在重構算法中加入全色參照信息提高光譜重構精度，但增加了系統復雜程度[21-22]。目前對重構算法的研究主要包括貪婪迭代算法、凸優化算法和基于貝葉斯框架的重構算法，前兩種算法在處理一維信號及二維圖像信號中取得了較好的成果，第三種算法由于其時間關聯性較強而不適用于光譜圖像的重構[23-26]。三維光譜數據具有維度高、矩陣規模大的特點。貪婪迭代算法在求解小規模一維和二維信號優化問題時獲得較好結果，但其在每一步都要完成信號或者殘差信號在冗余字典每一個原子上的投影計算，計算量成倍增加，重構時間較長，因而貪婪迭代算法不適用于混疊光譜數據重構求解。文獻[27]提出兩步迭代收縮閾值（Two-Step Iterative Shrinkage/Threshold，TwIST）算法實現CASSI系統光譜重構，但光譜重構精度較低。文獻[28]提出結合圖像自回歸模型的重構算法，因光譜圖像模型不同于二維圖像，該方法不適用于光譜數據重構。本文在分析TwIST算法基礎上，引入兩點改進：一是將求解函數中正則項設為光譜圖像的全變分正則約束項，使得求解結果最大可能的保留圖像本身的高頻信息；二是在每次迭代計算前，對更新項進行全變分去噪處理，避免引入的噪聲在后續迭代中放大。為了驗證基于全變分正則項的TwIST算法對CASSI數據具有較高的光譜重構精度，本文對仿真數據和試驗數據進行了重構。

1 CASSI數據構成

編碼孔徑快照式光譜儀數據可以表示為

編碼模板的空間調制使光譜數據損失部分空間信息，棱鏡色散后探測器同一像素上采集到不同空間位置光譜的混疊信息，因此重構過程是從不完整的觀測中恢復全部光譜數據，需求解嚴重的病態方程，可通過求解優化問題解決，目標函數表示為

式（2）右邊第一項為系統保真項，衡量優化結果與系統觀測之間的誤差；第二項為正則項，一般根據求解目標的固有性質對求解過程進行約束。

正則項對光譜重構精度有較大影響。稀疏重構的梯度投影算法以場景信息在空間上分段平滑作為先驗知識，實現優化算法的重構[29]?；谌兎值恼齽t約束項可以表示為

式中 （，）為離散后的空間坐標；為離散后的譜段；F,n,k，F+1,n,k，F,n+1,k分別為空間坐標為（,），（1,），（,+1）譜段的光譜數據。

本文采用全變分正則項作為求解優化問題的約束，在抑制噪聲的同時防止損失過多細節信息來保證光譜重構精度。

2 基于全變分正則項的TwIST算法

式（2）中的優化問題可由TwIST算法求解。本文通過改進TwIST算法，提高光譜重構精度。

2.1 TwIST算法

TwIST算法由迭代收縮閾值算法（Iterative Shrinkage/Thresholding，IST）和迭代加權收縮算法（Iterative Re-Weighted Shrinkage，IRS）結合而來，核心思想是利用前兩個估計值來更新當前值，即“兩步”迭代，主要迭代更新過程如下：

2.2 改進后的TwIST算法

全變分正則項基于光譜圖像空間分段平滑過渡的先驗知識[30]，對病態方程的求解施加約束。本文對TwIST算法提出兩點改進：一是將目標函數的約束項指定為光譜圖像的全變分正則約束項，使最終重構近似解保留圖像更多高頻細節，同時實現較快的收斂；二是在每次迭代計算前，對更新項進行全變分去噪處理，以免引入的噪聲隨迭代計算而放大。

式中（?）為光譜數據的目標函數。

圖1 改進的TwIST算法流程圖

3 仿真數據的重構及分析

3.1 建立仿真模型

CASSI系統數學模型為

式中 (,)表示連續的空間坐標；(,)為獲取的圖像；(,,)表示輸入的原始光譜信息；(,,)表示編碼模板空間調制和色散棱鏡光譜調制；(,)表示系統引入的各種噪聲。

離散形式可以表示為

通過Matlab實現CASSI成像的仿真建模。以隨機的0/1矩陣建立編碼模板的空間調制，采用理想線性色散模型。以高光譜數據集作為輸入，獲得的仿真圖像分辨率667像素×520像素，共30個譜段，譜段范圍426~670nm。圖2為仿真模擬采集結果。

圖2 仿真模擬采集結果

圖2（a）仿真目標場景包含黑白棋盤格、標準色卡以及其他目標物；圖2（b）仿真采集中光譜沿豎直方向色散，獲取含光譜信息的編碼混疊圖。

圖3為等效觀測矩陣仿真結果，共30個譜段，譜段范圍426~670nm。等效觀測矩陣在CASSI系統中起到空間調制作用，在算法中作為等效觀測矩陣參與光譜數據立方體重構。

圖3 等效觀測矩陣仿真結果

3.2 仿真光譜數據重構結果

圖4是以文獻[27]中的TwIST算法重構的仿真光譜數據。在426~496nm和587~670nm的光譜切片信息因噪聲等影響而損失嚴重，光譜數據整體的重構品質較差[13,27]。經計算平均光譜重構精度為73.21%。綜合不同仿真數據的重構結果，TwIST算法的光譜重構精度低于80%。

圖4 TwIST算法重構光譜數據

圖5為改進后的基于全變分正則約束的TwIST算法重構出的30個光譜切片。

圖5 基于全變分正則項的TwIST算法仿真數據重構光譜

對重構后的光譜切片進行分析可以看出，圖2（a）右側的色卡在圖5的光譜切片中清晰可見，色卡上不同塊的亮暗表示當前譜段下不同顏色灰度值的強度；圖2（a）左側的黑白色棋盤格在圖5不同光譜切片中邊緣清晰可見，重構的光譜數據較好地保留了空間細節信息，驗證了基于全變分正則項的TwIST算法對CASSI數據重構的可行性。

圖6為重構結果分析。從30個譜段中選取紅綠藍波段合成彩色圖，見圖6（a），與圖2（a）仿真目標場景相比，圖像的色彩和空間信息均具有較好的一致性。圖6（b）為場景中重構光譜曲線與原始光譜曲線的對比，重構光譜曲線與參照數據較好的吻合，驗證了基于全變分正則項的TwIST算法具有較好的光譜重構精度。

圖6 重構結果分析

選取圖2（a）中色卡作為分析目標，色卡有24個不同顏色，編號如圖7。

圖7 色卡的分塊編號

計算色卡各分塊的光譜重構精度

式中Q為重構數據譜段在空間位置的灰度值與總灰度值比值；F為重構數據譜段在空間位置的灰度值；D表示譜段在空間位置的灰度值與總灰度值比值；C表示譜段在空間位置的灰度值。

根據式（13）、（14）、（15）計算圖5各光譜切片色卡分塊的光譜重構精度，見表1。

表1 仿真色卡分塊重構精度

Tab.1 The reconstruction accuracy of simulation color blocks

根據表1可計算得出平均光譜重構精度為90.93%，優于TwIST算法的光譜重構精度，與輸入數據具有較好的一致性。

4 試驗數據的重構及分析

為了進一步驗證改進后的TwIST算法對CASSI數據的光譜重構精度，需要搭建CASSI試驗樣機，并采集數據。圖8為試驗采集結果，圖8（a）試驗目標場景，包含黑白棋盤格、標準色卡以及其他目標物，圖8（b）CASSI系統試驗采集圖，試驗采集中光譜沿豎直方向色散，展示了編碼后光譜的混疊結果。

圖8 試驗采集結果

圖9為等效觀測矩陣采集結果。圖9（a）為467nm譜段等效觀測矩陣的采集結果，圖9（b）為等效觀測矩陣局部細節。

圖9 等效觀測矩陣采集結果

圖10為基于全變分正則項的TwIST算法重構出的光譜數據切片，共有34個譜段，譜段范圍415~670nm。

圖10 基于全變分正則項的TwIST算法試驗數據重構光譜

對重構出的光譜切片進行分析可以看出，圖8（a）右側的色卡在圖10中不同光譜切片中清晰可見，色卡分塊的亮暗表示當前波長不同顏色灰度值的強度；圖8（a）左側的黑白色棋盤格在圖10中不同光譜切片中邊緣清晰可見，重構的光譜數據較好的保留了空間細節信息，這些驗證了基于全變分正則項的TwIST算法對CASSI實測數據重構的有效性。

按式（12）、（13）、（14）計算圖10各光譜切片中色卡分塊光譜重構精度，見表2。

表2 試驗色卡分塊光譜重構精度

Tab.2 The reconstruction accuracy of experimental color blocks

根據表2可計算得到平均光譜重構精度為86.56%，優于TwIST算法的光譜重構精度，與試驗實測數據具有較好的一致性。

5 結束語

本文基于全變分正則項，對TwIST算法提出兩點改進，并通過對仿真數據和試驗數據的重構，驗證了改進的基于全變分正則項的TwIST算法具有較好的光譜重構精度。改進的基于全變分正則項的TwIST算法對仿真數據光譜重構精度為90.93%，對試驗數據的光譜重構精度為86.56%，均優于TwIST算法的光譜重構精度。該優化設計的光譜重構算法可為CASSI光譜重構精度提升提供參考。

[1] DONOHO D L. Compressed Sensing[J]. IEEE Transactions on Information Theory, 2006, 52(4): 1289-1306.

[2] CANDES E J, WAKIN M B. An Introduction to Compressive Sampling[J]. IEEE Signal Processing Magazine, 2008, 25(2): 21-30.

[3] KUTYNIOK G. Compressed Sensing: Theory and Applications[J]. Corr, 2012, 52(4): 1289-1306.

[4] ELAD M. Optimized Projections for Compressed Sensing[J]. IEEE Transactions on Signal Processing, 2007, 55(12): 5695-5702.

[5] BLUMENSATH T, DAVIES M E. Iterative Hard Thresholding for Compressed Sensing[J]. Applied & Computational Harmonic Analysis, 2009, 27(3): 265-274.

[6] FIGUERIREDO M A T, NOWAK R D, WRIGH S J. Gradient Projection for Sparse Reconstruction: Application to Compressed Sensing and Other Inverse Problems[J]. IEEE Journal of Selected Topics in Signal Processing, 2008, 1(4): 586-597.

[7] SHI G M, LIU D H, GAO D H, et al. Advances in Theory and Application of Compressed Sensing[J]. Acta Electronica Sinica, 2009, 37(5): 1070-1081.

[8] 石光明, 劉丹華, 高大化, 等. 壓縮感知理論及其研究進展[J]. 電子學報, 2009, 37(5): 1070-1081. SHI Guangming, LIU Dianhua, GAO Dahua. Advances in Theory and Application of Compressed Sensing[J]. ACTA Electronica Sinica, 2009, 37(5): 1070-1081. (in Chinese)

[9] WAGADARIKAR A, JOHN R, WILLETT R, et al. Single Disperser Design for Coded Aperture Snapshot Spectral Imaging[J]. Applied Optics, 2008, 47(10): 44-51.

[10] WAGADARIKAR A A, PITSIANIS N P, SUN X, et al. Single-shot Compressive Spectral Imaging with A Dual-disperser Architecture[J]. Optics Express, 2007, 15(21): 14013-14027.

[11] ARCE G R, BRADY D J, CARIN L, et al. Compressive Coded Aperture Spectral Imaging: An Introduction[J]. IEEE Signal Processing Magazine, 2014, 31(1): 105-115.

[12] WAGADARIKAR A A, PITSIANIS N P, SUN X, et al. Video Rate Spectral Imaging Using A Coded Aperture Snapshot Spectral Imager[J]. Optics Express, 2009, 17(8): 6368-6388.

[13] WAGADARIKAR A A, PITSIANIS N P, SUN X, et al. Spectral Image Estimation for Coded Aperture Snapshot Spectral Imagers[J]. Proceedings of SPIE - The International Society for Optical Engineering, 2008, 7076(36): 6824-6833.

[14] LLULL P, LIAO X, YUAN X, et al. Coded Aperture Compressive Temporal Imaging[J]. Optics Express, 2013, 21(9): 10526-10545.

[15] ARGUELLO H, ARCE G R. Colored Coded Aperture Design by Concentration of Measure in Compressive Spectral Imaging[J]. IEEE Transactions on Image Processing, 2014, 23(4): 1896-1908.

[16] KITTLE D. Design and Fabrication of an Ultraviolet-visible Coded Aperture Snapshot Spectral Imager[J]. Optical Engineering, 2012, 51(7): 071403.

[17] WANG Lizhi, LI Huan, GAO Dahua, et al. A Flexible Design for Coded Aperture Snapshot Spectral Imager[C]// SPIE Vol.9273, Optoelectronic Imaging and Multimedia Technology III, International Society for Optics and Photonics, 2014.

[18] 石光明, 劉丹華, 高大化. 隨機編碼感知的高分辨率遙感光譜計算成像[J]. 航天返回與遙感, 2011, 32(5): 60-66. SHI Guangming, LIU Danhua, GAO Dahua. High-resolution Computational Spectral Imaging of Remote Sensing Based on Coded Sensing[J]. Spacecraft Recovery & Remote Sensing, 2011, 32(5): 60-66. ( in Chinese)

[19] 李國. 基于雙通道的壓縮光譜成像及其重構算法GPU實現[D]. 西安: 西安電子科技大學, 2014. LI Guo. The Double-channel Compressed Spectral Imaging and Reconstruction Algorithm Implemented on GPU[D]. Xi'an: Xidian University. 2014. (in Chinese)

[20] RUEDA H, ARGUELLO H, ARCE G R. Compressive Spectral Imaging Based on Colored Coded Apertures[C]// 2014 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), May 4-9, 2014, Florence, Italy. IEEE, 2014.

[21] WANG Lizhi，GAO Dahua，SHI Guangming，et al. Double-channel Compressive Spectral Imaging via Complementary Code Patterns[C]// 2013 IEEE International Conference on Signal Processing，Communication and Computing (ICSPCC), Aug 5-8, 2013, Kunming, China. IEEE, 2013: 1-5.

[22] WANG Lizhi, XIONG Zhiwei, GAO Dahua. et al. Dual-camera Design for Coded Aperture Snapshot Spectral Imaging[J]. Applied Optics, 2015, 54(4): 848-858.

[23] TROPP J A, GILBERT A C. Signal Recovery from Random Measurements Via Orthogonal Matching Pursuit[J]. IEEE Transactions on Information Theory, 2007, 53(12): 4655-4666.

[24] BOYD, VANDENBERGHE, FAYBUSOYICH. Convex Optimization[J]. IEEE Transactions on Automatic Control, 2006, 51(11): 1859-1859.

[25] CHAE B G, LEE S. Sparse-View CT Image Recovery Using Two-step Iterative Shrinkage-thresholding Algorithm[J]. ETRI Journal, 2015, 37(6): 1251-1258.

[26] LI Y, HU B L, WANG S, et al. Coded Aperture Spectral Data Restoration Based on Two-step Iterative Shrinkage/Threshold[J]. Spectroscopy and Spectral Analysis, 2014, 34 (3): 847-850.

[27] BIOUCAS-DIA J M, FIGUEIREDO M A T. A New TwIST: Two-step Iterative Shrinkage/Thresholding Algorithms for Image Restoration[J]. IEEE Transactions on Image Processing, 2007, 16(12): 2992-3004.

[28] WU Xiaolin, ZHANG Xiangjun, WANG Jia. Model-guided Adaptive Recovery of Compressive Sensing[C] //2009 Data Compression Conference (DCC 2009), March 16-18, 2009, Snowbird, UT, USA. IEEE, 2009: 123-132.

[29] XU H K. Averaged Mappings and the Gradient-projection Algorithm[J]. Journal of Optimization Theory and Applications, 2011, 150(2): 360-378.

[30] JIA R Q, ZHAO H. A Fast Algorithm for The Total Variation Model of Image Denoising[J]. Advances in Computational Mathematics, 2010, 33(2): 231-241.

A Reconstruction Algorithm of CASSI Data Based on Total Variation Regular Terms

WANG Yechao1,2CHEN Xiaoli1,2ZHONG Xiaoming1,2ZHAO Haibo1,2ZHANG Lisha1,2SU Yun1,2

（1 Beijing Institute of Space Mechanics & Electricity, Beijing 100094, China） （2 Key Laboratory for Advanced Optical Remote Sensing Technology of Beijing, Beijing 100094, China）

In designing spectral reconstruction algorithm for Coded Aperture Snapshot Spectral Imager (CASSI), the two-step iterative shrinkage/threshold (TwIST) algorithm can realize spectral reconstruction of CASSI, but the accuracy of spectral reconstruction is lower than 80% due to selection of regular terms and amplification of the introduced noise with the increase of iteration number. In this paper, based on the study of the traditional TwIST algorithm, two improvements are proposed, selecting the total variation regular constraint terms and denoising the updated terms in each iteration. In order to verify the improved algorithm, the data of CASSI are simulated by computer, with the spectral reconstruction accuracy 90.93%. The CASSI prototype is built to collect the experimental data, and the spectral reconstruction accuracy of the experimental data is 86.56%. The research results provide a reference for the following CASSI data reconstruction.

two-step iterative shrinkage/threshold (TwIST); compressed sensing; coded aperture; total variation; spectral reconstruction; spectral imager

TP751

A

1009-8518(2020)01-0091-11

10.3969/j.issn.1009-8518.2020.01.011

2019-10-16

北京市科技計劃課題（Z181100003018003）

王業超, 陳曉麗, 鐘曉明, 等. 基于全變分正則項的CASSI數據重構算法[J]. 航天返回與遙感, 2020, 41(1): 91-101.

WANG Yechao, CHEN Xiaoli, ZHONG Xiaoming, et al. A Reconstruction Algorithm of CASSI Data Based on Total Variation Regular Term[J]. Spacecraft Recovery & Remote Sensing, 2020, 41(1): 91-101. (in Chinese)

王業超，男，1993年生，2016年獲北京航空航天大學校飛行器設計專業學士學位，現在中國空間技術研究院飛行器設計專業攻讀碩士學位。研究方向為計算光譜成像技術。E-mail：cast_wangyc_508@163.com。

（編輯：王麗霞）