基于分塊對角矩陣的二維壓縮感知數據采集和重構方法

程 濤

(1.深圳大學光電子器件與系統重點實驗室，廣東 深圳 518060； 2.廣西科技大學汽車與交通學院，廣西 柳州 545006； 3.深圳生物醫學光學微納檢測與成像重點實驗室，廣東 深圳 518060)

基于分塊對角矩陣的二維壓縮感知數據采集和重構方法

程 濤1,2,3

(1.深圳大學光電子器件與系統重點實驗室，廣東 深圳 518060； 2.廣西科技大學汽車與交通學院，廣西 柳州 545006； 3.深圳生物醫學光學微納檢測與成像重點實驗室，廣東 深圳 518060)

針對遙感中線陣推掃數據采集模式和現有方法重構能力的不足，提出基于分塊對角矩陣和TV (Total Variation)算法的二維壓縮感知模型和方法。該方法能夠使壓縮感知的約束函數和目標函數同時包含完整的圖像二維信息，是真正完全意義上的二維壓縮感知。在不改變傳統數據采集模式的基礎上，通過分塊對角矩陣的后處理，實現二維壓縮感知。實驗結果表明，該方法使圖像的重構效果獲得極大改善，SNR提高約8 dB。但是該方法不適用于OMP(Orthogonal Matching Pursuit)和BP(Basis Pursuit)算法。該方法促進了定向遙感的發展，如與其他模型和方法結合，能進一步提高重構效果。

二維壓縮感知；定向遙感；線陣推掃；分塊對角矩陣；全變分

0 引言

根據奈奎斯特采樣定理， 采樣頻率是信號最高頻率2倍以上時，才能確保由采樣完全重構原始信號。壓縮感知能以遠低于奈奎斯特采樣定理要求的頻率采樣，在采集信號的同時實現對信號的壓縮，并能高質量地恢復原始信號[1-2]。因此，壓縮感知一出現就在圖像處理、視頻分析、雷達遙感、信息通信和醫學成像等領域成為研究熱點[3-4]。但是當前研究多是將壓縮感知作為一種數據后處理技術。真正能將壓縮感知應用于數據采集領域也只有CT(Computed Tomography，即電子計算機斷層掃描)和MRI(Magnetic Resonance Imaging，磁共振成像)等醫學成像和超分辨顯微成像等幾個能夠實現凝視狀態的領域[5-7]。

遙感數據的獲取一般分為星上和地面兩個階段。在衛星上需要采集、壓縮和傳輸數據，在地面需要接收、儲存、解壓和使用數據。遙感技術的發展使覆蓋全球海量遙感數據的獲取成為現實。如能實現采集壓縮一體化的信號數據采集處理格局，就可采集、傳輸、存儲、處理和管理很少的數據，從而擺脫傳統技術的窘境，節約巨大的人力物力資源。壓縮感知技術的出現為改變傳統的遙感數據獲取模式提供了可能。文獻[8-9]基于線陣推掃模式提出定向遙感和定向變化檢測。一般，采集2～3倍的測量數據就可完全重構變化區域。但是當前的各種壓縮感知重構算法和模型無法在線陣推掃模式下充分利用已有的先驗信息，因此重構效果并不好。本文針對當前遙感中基于線陣推掃模式的數據采集方法不能充分利用信號先驗信息的不足，提出了基于分塊對角矩陣和TV算法的二維壓縮感知模型和方法。

1 壓縮感知理論

壓縮感知理論表明，如果信號是稀疏或可壓縮的，就能以遠低于奈奎斯特(Nyquist)采樣定理的采樣率采集信號，并能以高概率精確或近似重建信號，壓縮感知模型如式(1)所示，求解式(1)就可以重構出信號x[10-11]。

min‖x‖0s.t.y=Φx

(1)

式中，y是測量數據，y∈RM；Φ是測量矩陣，Φ∈RM×N，M

測量矩陣Φ的性質是關系信號重構效果的關鍵。當前判斷測量矩陣Φ性能優劣的主要判據為是否滿足RIP(Restricted Isometry Property，約束等距性質)，如式(2)所示[10-11]。隨機矩陣的RIP性質一般都較好。

(1-σK)‖x‖2≤‖Φx‖2≤(1+σK)‖x‖2

(2)

式中，σK∈[0,1)。

如果測量矩陣滿足RIP和Johnson- Lindenstrauss定理，那么可導出RCP(Restricted Conformal Property，約束保角性質)，如式(3)所示[12]：

(3)

RIP主要反映了測量數據和信號之間的能量繼承關系，RCP反映了測量數據和信號之間的方向繼承關系。在二維情況下，這種繼承關系表現的非常明顯，并可作為判斷能量和方向的先驗信息[8-9,13]。

2 基于線陣推掃模式的二維壓縮感知

當前基于壓縮感知的遙感研究多是基于一維信號。即使研究圖像這種二維信號也多是轉化為一維信號后再作壓縮感知研究。在壓縮感知的各種重構算法中，只有最小全變分法(Total Variation，TV)在目標函數中利用了圖像像素間的灰度梯度信息，但在約束函數中依然把二維信號一維化。從重構效果看，只有在梯度大的地方優于離散余弦變換(Discrete Cosine Transform, DCT)?；叶忍荻染褪沁b感影像的結構之一，是很好的先驗信息，如充分利用可有助于改善提高信號的重構精度和效率。但把二維圖像一維化破壞了這種結構先驗信息[8,14]。

遙感衛星等的影獲取多是采用線陣推掃模式，設掃描條帶的長度為L。如果以矩陣X(X∈RN×L)表示掃描條帶，則衛星或飛機的掃描線每次只能取得X某列的數據。以二維形式表示其壓縮感知模型，則如式(4)所示：

min‖xj‖0s.t.Y=ΦX,j∈[1,L]

(4)

式中，Y∈RM×L，xj是X的列向量。由式(4)可見，二維化后的測量矩陣Φ的規模為M×N；如果將遙感影像轉化為一維信號，測量矩陣Φ的規模為ML×NL，前者僅為后者的1/L2。測量矩陣規模的縮小將會大大提高數據解算的速度和效率。式(4)也可等價表示成式(5)，兩種表示方法并無本質區別，因此式(4)并不是完全意義上的二維壓縮感知模型。但是式(4)的約束函數能更好地反映信號和測量數據的二維信息[8]。

min‖xj‖0s.t.yj=Φxj

(5)

文獻[13]基于式(5)采用OMP算法逐列重構，然后通過RCP，利用二維圖像的列間相關性，對重構結果修正，取得了更好的重構效果。對于可壓縮信號，TV是效果最好的重構算法。文獻[9]對滿足稀疏要求的列采用OMP算法重構，然后利用RCP先驗信息確定出可壓縮列，并對之采用TV重構，取得了很好的重構效果。這些方法在壓縮感知模型的目標函數中并沒有引入列間的相關性信息，而是通過后處理的形式判斷重構結果的優劣，并對重構不好的區域修正。因此本質上都不是真正的二維壓縮感知。

min‖X‖TVs.t.Y=ΦX

(6)

如能將式(4)轉化為式(6)的形式就能實現真正的二維壓縮感知。因為式(6)的目標函數既包含了圖像縱向(列向)的灰度梯度信息，也包含了圖像橫向的灰度梯度信息。約束函數也完整保留了圖像和測量數據的二維結構信息。

但是式(6)并不能直接應用于現有壓縮感知的TV模型。因此，如將式(6)中的圖像，測量數據都一維化，那么式(6)就等價轉化為式(7)。式(6)和式(7)的目標函數盡管寫法不一樣，但本質上是完全等價的，都包含了圖像橫向和縱向的全部梯度信息；約束函數在數學上也是完全等價的，并沒有導致信息的丟失和減少。

min‖x‖TVs.t.y=ΦDIGx

(7)

3 實驗結果分析和驗證

所有的實驗采用Matlab在普通臺式機上完成，CPU主頻3.1 GHz，4核，內存1 GB。Matlab的操作界面如圖2所示。用灰度圖表示的子塊矩陣Φ、分塊對角矩陣ΦDIG和稠密矩陣ΦDEN如圖3所示。子塊矩陣Φ的矩陣規模為M×N；分塊對角矩陣ΦDIG的矩陣規模為ML×NL，其中包含L個子塊矩陣Φ；稠密矩陣ΦDEN的矩陣規模為ML×NL，其大小是子塊矩陣Φ的L2倍。子塊矩陣Φ和稠密矩陣ΦDEN都是服從N(0，1/M)分布的高斯測量矩陣，兩者的區別只在于矩陣大小規模的不同。

基于子塊矩陣Φ的重構算法流程，就是通過Φxj逐列采集測量數據yj；然后，基于式(5)，通過TV算法逐列重構xj；最后，將逐列重構得到的所有xj，排列組合成二維圖像X。

基于分塊對角矩陣ΦDIG的重構算法流程，就是通過如圖1所示的算法流程圖重構得到二維圖像X。

基于稠密矩陣ΦDEN的重構算法流程，就是通過稠密矩陣ΦDEN直接采集二維圖像X的一維化數據x，ΦDENx；然后，基于min‖x‖0s.t.y=ΦDENx，通過TV算法重構得到x；最后將x二維化成二維圖像X。

圖4和表1是基于子塊矩陣Φ、分塊對角矩陣ΦDIG和稠密矩陣ΦDEN對Boats、Cameraman、Lena和遙感影像Mulargia lake以及墨西哥森林5幅圖像采用TV算法的重構結果。Mulargia lake是意大利撒丁島Mulargia湖的兩個時相遙感影像的差值圖，由Landsat-5衛星在波段4分別拍攝于1995年09月和1996年07月。該圖像反映了Mulargia湖水位上升而造成的陸地淹沒情況。為滿足定向遙感的稀疏條件，對兩個時相的原始遙感影像做簡單處理，使其差值圖中未變化區域的灰度值都為0。墨西哥森林是拍攝于某時刻的遙感圖像，背景與目標灰度非常相近，幾乎無法區分目標和背景。圖4中，除Mulargia Lake外的其他4幅圖都是可壓縮的。Mulargia Lake圖中的黑色區域的像素灰度都為0。除此之外的其他區域也是可壓縮的。表1中，包含橫向灰度梯度信息的目標函數值用TV表示，橫向灰度梯度信息用TVx表示，TV和TVx的值×104后才是其真實值。

子塊矩陣基于式(5)逐列重構，分塊對角矩陣基于式(7)重構。稠密矩陣也基于式(7)重構，不過這時式(7)中的約束函數的分塊對角矩陣用稠密矩陣替換。子塊矩陣和稠密矩陣都是高斯矩陣。

圖4中自左向右依次是基于子塊矩陣、分塊對角矩陣和稠密矩陣的重構結果，最右側是原始圖像。對比左側的兩列圖像可以發現，基于TV算法的分塊對角矩陣的重構效果遠遠好于子塊矩陣逐列重構的效果。由表1可見，前4幅圖的SNR提高約8 dB，計算時間最大約增加22倍(通過Matlab中的tic和toc命令統計計算時間)。而對于背景與目標灰度非常相近的墨西哥森林，SNR僅提高約2.5 dB，計算時間約增加14倍。包含橫向灰度梯度信息的最終目標函數數值(TV)和橫向灰度梯度值(TVx)都小于子塊矩陣，更接近于稠密矩陣。

SNR計算公式如式(8)所示：

(8)

式中，x是真實信號，x∈RN；‖·‖2表示向量的模；xR表示x的重構結果。

對比左側的兩幅Lena圖可以發現，Lena圖中左側立柱的重構效果，分塊對角矩陣劣于子塊矩陣。這說明，逐列重構的子塊矩陣更適用于灰度梯度縱向變化平緩的圖像局部區域。對比左側的兩幅Cameraman圖可以發現，在頭部附近的縱向區域，逐列重構的子塊矩陣的重構效果也劣于分塊對角矩，這是因為該區域的縱向灰度梯度變化大。分塊對角矩陣更適用于灰度梯度橫向變化平緩的圖像局部區域；例如，Boats圖的橫向灰度梯度變化較小，縱向灰度梯度變化大，所以子塊矩陣逐列重構得到的左側Boats圖的所有區域都很模糊。

表1 3類矩陣的SNR和計算時間

如果將逐列重構的子塊矩陣和分塊對角矩陣相結合，在縱向灰度梯度變化平緩的區域采用逐列重構的子塊矩陣，在其他區域使用分塊對角矩陣就能取得更好的重構效果。

由圖4和表1可見，稠密矩陣的重構效果和SNR盡管都好于分塊對角矩陣，但是稠密矩陣這種矩陣形式并不能應用于基于線陣推掃的遙感數據采集模式中。由表1可見，分塊對角矩陣和稠密矩陣的計算速度幾乎相當，遠遠慢于子塊矩陣。這是因為子塊矩陣的矩陣規模遠遠小于分塊對角矩陣和稠密矩陣。

由于在TV算法中， 3類矩陣的目標函數包含信息量的不同，導致無法準確測試3類矩陣的性能和效果，因此，基于OMP和BP算法測試3類矩陣的性能?；诓煌∈瓒鹊母咚剐盘枺鲗嶒?4次，計算精確重構概率。圖5和表2是基于子塊矩陣、分塊對角矩陣和稠密矩陣對同一套測試數據分別采用OMP和BP算法的重構結果。由圖5中OMP和BP可見，稠密矩陣在稀疏度較小時重構效果好于子塊矩陣；當稀疏度較大時反之；分塊對角矩陣的重構效果最差。

由表2可見，OMP算法，子塊矩陣的計算速度最快。分塊對角矩陣和稠密矩陣的計算速度相近，都遠遠慢于子塊矩陣。BP算法，稠密矩陣的計算速度最慢。子塊矩陣和分塊對角矩陣的計算速度相近，都遠遠快于稠密矩陣。

時間/s子塊矩陣分塊對角矩陣稠密矩陣OMP算法0 090215 382224 746BP算法4 4693 5644236 520

在TV算法中，分塊對角矩陣和稠密矩陣的目標函數不但包含原始圖像的縱向灰度梯度信息，而且還包含橫向的灰度梯度信息；但是子塊矩陣只包含原始圖像的縱向灰度梯度信息。OMP和BP算法表明，分塊對角矩陣性能劣于子塊矩陣和稠密矩陣。但是在TV算法中，包含橫向灰度梯度信息的目標函數彌補，并遠遠超過了分塊對角矩陣性能的不足。

文獻[9]說明TV算法只適用于可壓縮信號，對于稀疏信號，使用其他算法可取的更好的重構效果。因此，在定向遙感的實際使用中，對稀疏區域使用其他重構算法重構，對于可壓縮區域中縱向灰度梯度變化平緩的使用逐列重構的子塊矩陣；對于可壓縮區域中縱向灰度梯度變化不平緩區域的使用分塊對角矩陣重構。

4 結論

本文提出了基于分塊對角矩陣和TV算法的二維壓縮感知模型和方法。該方法能夠使壓縮感知的約束函數和目標函數同時包含完整的圖像二維信息，是真正完全意義上的二維壓縮感知。實驗結果表明，由于目標函數中包含了圖像的二維梯度信息，該方法使圖像的重構效果獲得極大改善，SNR提高約8 dB，但是該方法不適用于OMP和BP算法。這是因為它們的目標函數中并不包含圖像的二維梯度信息，本質上依然是一維壓縮感知。利用更多圖像像素相關信息的目標函數能改善重構效果。該方法促進了定向遙感的發展，擬結合其他模型和方法，進一步提高重構效果。

[1]CANDES E J, ROMBERG J, TAO T. Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information [J]. Information Theory, IEEE Transactions on, 2006, 52(2): 489-509.

[2]CANDES E J, TAO T. Decoding by linear programming [J]. Information Theory, IEEE Transactions on, 2005, 51(12): 4203-4215.

[3]焦李成, 楊淑媛, 劉芳, 等. 壓縮感知回顧與展望 [J]. 電子學報, 2011, 39(7): 1651-1662.

[4]林杰, 石光明, 董偉生. 基于信息自由度采樣的信號重構方法研究進展 [J]. 電子學報, 2012, 40(8): 1640-1649.

[5]ZHU L, ZHANG W, ELNATAN D, et al. Faster STORM using compressed sensing [J]. Nature methods, 2012, 9(7): 721-723.

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

[8]程濤, 朱國賓, 劉玉安. 基于二維壓縮感知的定向遙感和變化檢測研究 [J]. 紅外與毫米波學報, 2013, 32(5): 491-496.

[9]CHENG T. Directional Remote Sensing [J]. Geod list, 2015, 52(4): 251-262.

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

[11]CANDES E J. Compressive sampling [C]//Proceedings of the International Congress of Mathematicians. Madrid, Spain: European Mathematical Society, 2006: 1-20.

[12]CHENG T. Restricted Conformal Property of Compressive Sensing [C]//2014 11th International Conference on Wavelet Active Media Technology and Information Processing. Chengdu: UESTC Press, 2014: 152-161.

[13]程濤, 朱國賓, 李小龍. 基于信號方向信息的二維壓縮感知影像重構研究 [J]. 半導體光電, 2014, 35(4): 32-37.

[14]程濤. 基于壓縮感知的遙感變化檢測研究 [D]. 武漢: 武漢大學, 2013.

Two-dimensional Compressive Sensing Data Acquisition and Reconstruction Based on Block Diagonal Matrix

CHENG Tao1,2,3

(1.Key Laboratory of Optoelectronic Devices and Systems, Shenzhen University, Shenzhen 518060, China; 2.Automotive & Transportation Engineering Institute，Guangxi University of Science and Technology, Liuzhou 545006, China; 3.Shenzhen Key Laboratory of Micro-Nano Measuring and Imaging in Biomedical Optics, Shenzhen University, Shenzhen, 518060, China)

For the data acquisition modes of linear array push-broom in remote sensing and deficiencies of existing methods for data reconstruction, the two-dimensional compressive sensing models and methods based on a block diagonal matrix and TV (Total Variation) algorithm were proposed in this paper. This method enabled the constraint and objective functions of compressive sensing to contain the completed two-dimensional image information. And it was totally real two-dimensional compressive sensing. On the basis of the traditional data acquisition mode, the two dimensional compressive sensing was realized by the post processing of the block diagonal matrix. Experimental results showed that the image reconstruction effect was improved greatly, and the SNR increased by about 8dB. However, this method was not suitable for OMP (Orthogonal Matching Pursuit) and BP (Basis Pursuit).

two-dimensional compressive sensing; directional remote sensing; linear array push-broom; block diagonal matrix; TV(total variation)

2017-01-10

國家自然科學基金項目資助(41461082,81660296); 中國博士后科學基金項目資助(2016M592525);廣西自然科學基金項目資助(2014GXNSFAA118285);廣西高校科學技術研究項目資助(YB2014212);廣西科技大學博士基金項目資助(?？撇?3Z12)

程濤(1976—)，男，廣西柳州人，博士，副教授，研究方向：壓縮感知和遙感研究。E-mail:ctnp@163.com。

TP391

A

1008-1194(2017)03-0060-06