解成俊，徐 林，張鐵山

（北華大學(xué) 數(shù)字圖像處理研究所，吉林 吉林 132011）

壓縮感知[1]（Compressed Sensing，CS）理論從某種意義上突破了香農(nóng)采樣定理和奈奎斯特采樣頻率的限制，CS理論帶來了信號采樣理論的變革，具有廣闊的應(yīng)用前景[2]。利用CS技術(shù)采集圖像數(shù)據(jù)可以大大降低數(shù)碼率[3]。美國RICE大學(xué)的Baraniuk等人運(yùn)用CS理論成功研制了單像素相機(jī)[4]。由于該相機(jī)直接獲取的是M次隨機(jī)線性測量值，而不是獲取原始信號的N（M?N）個像素值，因此為低碼率相機(jī)拍攝高質(zhì)量圖像提供了可能。CS理論也被應(yīng)用于雷達(dá)成像領(lǐng)域，與傳統(tǒng)雷達(dá)成像技術(shù)相比壓縮感知雷達(dá)成像實現(xiàn)了兩個重要改進(jìn)：在接收端省去脈沖壓縮匹配濾波器，同時由于避開了對原始信號的直接采樣，降低了接收端對模數(shù)轉(zhuǎn)換器件帶寬的要求，設(shè)計重點(diǎn)由傳統(tǒng)的設(shè)計昂貴的接收端硬件轉(zhuǎn)化為設(shè)計新穎的信號重構(gòu)算法，從而簡化了雷達(dá)成像系統(tǒng)[5]。Bhattacharya等還將CS理論應(yīng)用到合成孔徑雷達(dá)圖像數(shù)據(jù)獲取上[6]，解決了海量數(shù)據(jù)采集和存儲問題，顯著降低了衛(wèi)星圖像處理的計算代價。壓縮感知理論也被應(yīng)用于醫(yī)學(xué)成像領(lǐng)域，如稀疏核磁共振成像[7]、壓縮傳感三維磁共振波譜成像[8]。壓縮感知理論研究的核心問題是隨機(jī)測量和重構(gòu)問題。在信號重構(gòu)算法中基追蹤算法[9]是基本算法，該算法是從完備的(過完備)函數(shù)集合中得到信號的最稀疏表示，即用盡可能少的基精確地表示測量信號，要在測量矩陣的所有列向量中極小化一個全局目標(biāo)函數(shù)，其計算量仍然是巨大的[10]，因此基追蹤算法基本用在一維信號處理，為了用基追蹤算法實現(xiàn)求解二維問題，需對該算法進(jìn)行改進(jìn)，即修改基追蹤算法的約束條件問題變?yōu)榍蠼馄鋬?yōu)化問題,同時進(jìn)行了實驗測試，證明本文算法是有效的。

1 壓縮感知理論

壓縮感知理論是一種全新的在采樣的同時實現(xiàn)壓縮目的的理論框架，壓縮感知主要包括3個環(huán)節(jié)：信號稀疏化、隨機(jī)測量和重構(gòu)。設(shè)原始信號z∈RN，N是信號長度，Ψ是原始信號z的稀疏表示，x=Ψz，信號隨機(jī)測量的目的就是得到M（M?N）個觀測值 y=Φx=ΦΨz，y∈RM，其中 x表示長度為 N的k-稀疏信號 （K?N），k被稱為稀疏度，Φ∈RM×N為隨機(jī)測量矩陣。由于M?N，使得y=Φx=ΦΨz是欠定的，存在多個解，可以用式（1）求解 x。

但是問題（1）是一個NP難題，因此通常將其轉(zhuǎn)化為等效的l1凸集優(yōu)化問題[9]，即：

本文用改進(jìn)的基追蹤算法對式（2）進(jìn)行求解，再通過稀疏化的逆過程求出原始信號z。

2 改進(jìn)的基追蹤算法

由于測量矩陣Φ的設(shè)計滿足受限等距映射特性(Restricted Isometry Property，RIP)[11]，則重建信號x?的誤差可以表示為[9]：

修改式（2）使其變?yōu)榍蠼馐剑?）的優(yōu)化問題。

其中，c1為正常數(shù)，R是范數(shù)變化的半徑，p控制收斂速度，p越小收斂速度越快。

根據(jù)文獻(xiàn)[10]，有 k≤cM/log（N/M），從而得

其中，r=1/p-1，由此可知在1-范數(shù)的情況下，誤差可以控制在一定的范圍內(nèi)。

為求解（4）式，令 x=a-b，a≥0，b≥0，且式中 ai=（xi）+，bi=（-xi）+，i=1，2，…，n，n 為向量的維數(shù)，在此（·）+是對數(shù)據(jù)‘·’的實部進(jìn)行運(yùn)算，令 l=[1，1，…，1]，l的維數(shù)為則x的1-范數(shù)可以寫成‖x‖1=lTu。u的維數(shù)為n。

同理可以令 y-Ax=e-f，e≥0，f≥0，且式中 ei=（（y-Φx）i）+，fi=（-（y-Φx）i）+，i=1，2，…，m，m 為向量的維數(shù)，令 l=[1，1，…，1]，l的維數(shù)為則y-Φx的1-范數(shù)可以寫成‖y-Φx‖1=lTv。v的維數(shù)為m。

因此可以寫出下式：

這里cTX是目標(biāo)函數(shù)，DX=y是約束條件，X≥0是邊界條件。而對式（5）的求解見文獻(xiàn)[12]。

3 實驗結(jié)果及討論

本文采用的測試圖像是國際標(biāo)準(zhǔn)測試圖像lena.bmp(256×256×8 bits),利用小波變換實施稀疏化，設(shè)小波變換的尺度用J表示，對J=1、2、3的情況分別進(jìn)行采樣，用f代表采樣率。利用標(biāo)準(zhǔn)偽隨機(jī)數(shù)均勻分布序列矩陣P、Q和中心傅里葉變換F生成隨機(jī)測量矩陣Φ=PFQ，通過Φ對稀疏化后的數(shù)據(jù)進(jìn)行隨機(jī)測量，利用改進(jìn)的BP算法對測量后的數(shù)據(jù)進(jìn)行圖像重建。表1是對測試圖像進(jìn)行1層小波變換并對子帶LH1、HL1、HH1按不同采樣率進(jìn)行采樣后進(jìn)行重建的實驗結(jié)果。表 2是固定了子帶 LH1、HL1、HH1的采樣率為 0.1%,即采樣子帶 LH1、HL1、HH1的 49 個 數(shù) 據(jù) 點(diǎn)， 變 化 子 帶 LH2、HL2、HH2的采樣率，利用改進(jìn)的BP算法進(jìn)行重建的實驗結(jié)果。表3 是固 定了子帶 LH1、HL1、HH1和 LH2、HL2、HH2的采樣 率分別為0.1%和30%,即采樣子帶LH1、HL1、HH1的49個數(shù)據(jù)點(diǎn)和子帶 LH2、HL2、HH2的 3 686 個點(diǎn)，改變子帶 LH3、HL3、HH3的采樣率，利用改進(jìn)的BP算法進(jìn)行重建的實驗結(jié)果。隨著J的增加，重建圖像的PSNR降低，其主要原因是小波逆變換傳遞誤差造成的，而BP算法重建圖像過程中又再次放大了傳遞誤差。

表1 J=1時對高頻子帶采樣重建圖像效果Tab.1 Reconstruction result of the high frequency subband under sampling when J=1

圖1是仿真實驗采用的標(biāo)準(zhǔn)測試圖像lena.bmp(256×256×8 bits)； 圖 2是一層小波變換后對高頻子帶 LH1、HL1、HH1采樣0.1%的數(shù)據(jù)，總采樣率ftotal=25.08%時重建圖像的效果；圖3是二層小波變換后對高頻子帶LH1、HL1、HH1采樣0.1%的數(shù)據(jù)，LH2、HL2、HH2采樣 50%的數(shù)據(jù)， 總采樣率 ftotal=15.63%時重建圖像的效果，圖4是三層小波變換后對高頻子帶 LH1、HL1、HH1采 樣 0.1%的 數(shù)據(jù) 、LH2、HL2、HH2采 樣 30%的數(shù)據(jù)，LH3、HL3、HH2采樣 5%的數(shù)據(jù)，總采樣率 ftotal=7.5%時重建圖像的效果。由實驗結(jié)果可知，采樣很少的數(shù)據(jù)通過壓縮感知重建仍能得到效果很好的重建圖像。

表2 J=2時對高頻子帶采樣重建圖像效果Tab.2 Reconstruction result of the high frequency subband under sampling when J=2

表3 J=3時對高頻子帶采樣重建圖像效果Tab.3 Reconstruction result of the high frequency subband under sampling when J=3

圖1 原始圖像Fig.1 Original image

圖2 J=1，采樣率：f（LH1、HL1、HH1）=0.1%ftotal=25.08%，重建圖像PSNR=31.16Fig.2 J=1，sampling rate：f（LH1、HL1、HH1）=0.1%ftotal=25.08%，PSNR of reconstruction image=31.16

圖3 J=2，采樣率：f（LH1、HL1、HH1）=0.1%，f（LH2、HL2、HH2）=50%，ftotal=15.63%，重建圖像PSNR=24.13Fig.3 J=2，sampling rate：f（LH1、HL1、HH1）=0.1%，f（LH2、HL2、HH2）=50%，ftotal=15.63%，PSNR of reconstruction image=24.13

4 結(jié) 論

采用改進(jìn)的基追蹤算法求解二維壓縮感知重建問題可以得到令人滿意的結(jié)果，尤其是低碼率壓縮感知重建具有更好的效果，為我們研制開發(fā)基于壓縮感知的新一代設(shè)備提供了理論支持。

圖4 J=3，采樣率：f（LH1、HL1、HH1）=0.1%，f（LH2、HL2、HH2）=30%，f（LH3、HL3、HH3）=5%，ftotal=7.5%，重建圖像 PSNR=22.81Fig.4 J=3，sampling rate：f（LH1、HL1、HH1）=0.1%，f（LH2、HL2、HH2）=50%，f（LH3、HL3、HH3）=5%， ftotal=7.5%，PSNR of reconstruction image=22.81

[1]Donoho D.Compressed sensing[J].IEEE Transactions on Information Theory，2006，52（4）：1289-1306.

[2]Candès E， Romberg J，Tao T.Robust uncertainty principles：Exact signal reconstruction from highly incomplete frequency information [J].IEEE Transactions on Information Theory，2006，52（2）：489-509.

[3]Romberg J.Imaging via compressive sampling[J].IEEE Signal Processing Magazine，2008，25（2）：14-20.

[4]Gemmeke J F，Cranen B.Using sparse representations for missing data imputation in noise robust speech recognition[C]//European Signal Processing Conf.（EUSIPCO），Lausanne，2008：987-991.

[5]Duarte M F，Davenport M A，Takhar D，et al.Single-Pixel imaging via compressive sampling[J].IEEE Signal Processing Magazine，2008，25（2）：83-91.

[6]Baraniuk R，Steeghs P.Compressive radar imaging[C]//Proceedings of the Radar Conference.Washington D.C.，2007 ，128-133.

[7]Bhattacharya S，Blumensath T，Mulagrew B，et al.Fast encoding of synthetic aperture radar raw data using compressed sensing[C]//In：Proceedings of Statistical Signal Processing.Washington D.C.，2007，448-452.

[8]Lustig M，Donoho D L，Pauly J M.Sparse MRI：The application of compressed sensing for rapid MR imaging[J].Magnetic Resonance in Medicine，2007，58（6）:1182-1195.

[9]Chen S S，Donoho D L，Saunders M A.Atomic decomposition by basis pursuit[J].SIAM Review，2001，43（1）：129-159.

[10]Candès E.Compressive sampling[C]//Proceedings of the Internation Congress of Mathematicians，Madrid，Spain，2006，：1433-1452.

[11]Candes E，Wakin M B.An introduction to compressive sampling[J].IEEE Signal P rocessing Magazine，2008，25（3）：21-30.

[12]Bloomfield P，Steiger W.Least absolute deviation：theory，applications and algorithms[M].Boston：Birkhauser，1983.