閆 歌，許廷發(fā)，馬 旭，張宇寒，王 茜，譚翠媚

(北京理工大學 光電學院 光電成像技術(shù)與系統(tǒng)教育部重點實驗室，北京 100081)

1 引 言

高光譜遙感成像可同時獲取高空間分辨率和高光譜分辨率的圖像，在礦物勘探、環(huán)境監(jiān)測和食品檢測等方面有著廣泛的應(yīng)用[1-6]。由于應(yīng)用面廣，高光譜成像設(shè)備的發(fā)展十分迅速，逐漸朝著輕便小型化的方向發(fā)展[7-8]。但高光譜圖像空間分辨率、光譜分辨率的不斷提高，導致圖像數(shù)據(jù)量成倍增加，給數(shù)據(jù)的存儲和傳輸增加了壓力，限制了高光譜成像技術(shù)的進一步發(fā)展。因此，如何有效地壓縮高光譜圖像是當前該領(lǐng)域的研究熱點。

壓縮感知(Compressive Sensing,CS)理論指出，如果信號在某一變換域是 K-稀疏的或者說是可壓縮的，那么可以設(shè)計一個與相變換基不相關(guān)的測量矩陣來對信號進行觀測。其可利用低維測量值，通過求解一個凸的最優(yōu)化問題來實現(xiàn)原始信號的重構(gòu)[9-10]。CS邊壓縮邊采樣的特性，相對于傳統(tǒng)數(shù)據(jù)傳輸省去了先獲得大量冗余數(shù)據(jù)再壓縮的麻煩，降低了采樣成本，減少采樣時間。CS主要分為信號壓縮感知測量和信號重構(gòu)兩個步驟，在壓縮過程中運算量較少，而重構(gòu)則需要經(jīng)過大量運算多次迭代獲得最優(yōu)解。因此大多數(shù)壓縮感知算法的研究重點都在重構(gòu)算法方面[11-13]。常見的重構(gòu)算法有凸優(yōu)化算法、貪婪算法等，其中凸優(yōu)化算法重構(gòu)精度高但計算量大，貪婪算法計算量小但重構(gòu)精度低。影響信號恢復的質(zhì)量除了重構(gòu)算法外還有其他因素，例如能否找到合適的稀疏基對信號進行稀疏表示、觀測矩陣能否在采集過程中保留大部分原始信號的信息，以及能否在壓縮測量過程中合理地分配采樣率等。在視頻壓縮感知(Compressive Video Sensing)中，由于圖像規(guī)模較大，需要先進行分塊，然后分別對每一個圖像塊進行壓縮測量。針對這一問題，已經(jīng)有不少學者提出了自適應(yīng)采樣的方案，例如文獻[14]利用圖像塊方差作為衡量細節(jié)復雜度的參數(shù)，為各塊設(shè)定不同的測量率，但由于不同的圖像方差塊差距可能巨大，導致分配的測量值差距也很大，甚至有可能出現(xiàn)采樣值超出圖像塊像素數(shù)的情況；Hung-Wei Chen等人[15]利用前一幀圖像作為參考幀，訓練出圖像的稀疏字典，通過訓練出的稀疏字典對參考幀進行稀疏表示，再根據(jù)圖像塊稀疏表示后的稀疏系數(shù)的方差自適應(yīng)分配采樣率。但假如圖像使用普通稀疏基(如DCT基，小波基等)來稀疏表示，不訓練稀疏字典的情況下，得到的各圖像稀疏系數(shù)方差的差距很大，并不能合理分配測量值；左覓文等人[16]提出了一種基于幀間相關(guān)性的自適應(yīng)采樣率分配方案，將圖像塊采樣分為固定預采樣、基于圖像塊復雜度比例的自適應(yīng)采樣兩個采樣部分，根據(jù)圖像塊的變化程度動態(tài)調(diào)整測量率。上述方案都是基于視頻信號的幀間相關(guān)性很高這一特性來進行自適應(yīng)采樣的。高光譜圖像相鄰波段間的相關(guān)性也很高，目前的利用高光譜圖像光譜相關(guān)性的研究主要在譜間預測重構(gòu)研究方面[17-18]，并且目前的高光譜壓縮感知算法普遍以固定的測量率對每個波段的圖像塊進行測量[19]。用高光譜圖像強譜間相關(guān)性，可以將自適應(yīng)采樣率應(yīng)用在高光譜圖像壓縮感知領(lǐng)域。另外，由于目前高光譜壓縮感知成像系統(tǒng)只能同時對所有圖像塊采集相同測量值，上述的自適應(yīng)測量方法雖然提高了圖像的重構(gòu)精度，但是由于存在方差懸殊等問題，導致不同圖像塊分配的測量值數(shù)目也相差很大，這對自適應(yīng)壓縮感知算法應(yīng)用于高光譜成像系統(tǒng)帶來了困難。

本文提出一種利用邊緣信息動態(tài)測量高光譜圖像的壓縮感知算法。該方法有效利用高光譜圖像波段間的相關(guān)性，根據(jù)各波段圖像塊含有邊緣信息的多少自適應(yīng)地調(diào)整測量率，在給定總采樣率的條件下可實現(xiàn)圖像塊的采樣率動態(tài)分配。由于在背景區(qū)域?qū)γ總€圖像塊都給定了最少的測量值數(shù)目，因此既不會犧牲背景的重構(gòu)效果，又控制了圖像塊間測量值數(shù)目差異，使自適應(yīng)壓縮感知理論應(yīng)用在系統(tǒng)中成為可能。

2 動態(tài)測量的高光譜圖像壓縮感知

為了獲得邊緣信息，就要先獲取原始信號?；诟吖庾V圖像譜間相關(guān)性很高的特性，可以先用傳統(tǒng)壓縮感知方法用固定測量率的隨機矩陣對某個波段先進行觀測，獲得測量值后重構(gòu)出該波段圖像，作為其他波段的先驗信息。本文提出的利用邊緣信息的動態(tài)測量壓縮感知高光譜圖像重構(gòu)算法框架，如圖1所示。在該框架中，對參考波段采用傳統(tǒng)的固定采樣率測量，然后重構(gòu)出參考波段的圖像，獲得圖像的先驗信息后，再對非參考波段圖像采用本文提出的動態(tài)測量率進行采樣與壓縮感知重構(gòu)。

在對參考波段圖像進行測量時，每個圖像塊的測量率是固定的，通過傳統(tǒng)壓縮感知方法測量并重構(gòu)出圖像后，對該圖像進行邊緣提取，獲得邊緣信息。由于高光譜圖像譜間相關(guān)性很高，并且圖像在空間方向上是不移動的，那么提取的邊緣區(qū)域就可以看作是所有波段的邊緣區(qū)域。有了高光譜圖像的邊緣信息后，就可以根據(jù)待測量的圖像塊中含有邊緣信息的多少來分配測量值，含有邊緣信息多的圖像塊分配的采樣率高，含有邊緣信息少的圖像塊分配的采樣率低。

圖1 利用邊緣信息的動態(tài)測量壓縮感知高光譜圖像重構(gòu)算法 Fig.1 CS hyperspectral image reconstruction algorithm based on edge information

2.1 高光譜圖像參考波段的觀測

對高光譜圖像參考波段的觀測是采用固定測量率的分塊壓縮感知方法。參考波段Nr×Nc個像素的圖像Xkey包含了n個大小為B×B的不重疊塊，第i個圖像塊的列向量表示為xi,i=1,2,…,n，通過對每個圖像塊用相同的隨機投影矩陣觀測，可以得到第i個圖像塊的壓縮測量值

yi=Φxi,

(1)

式中，Φ是一個大小為M×B2的測量矩陣，M為測量次數(shù)，并有M<

R=M/B2,

(2)

參考波段的重構(gòu)圖像將用來提取邊緣區(qū)域，然后利用邊緣信息作為分配動態(tài)測量率的重要參數(shù)。參考波段的重構(gòu)質(zhì)量將影響邊緣提取的準確性，從而影響整個高光譜圖像的重構(gòu)效果。所以，參考波段會采用較高的采樣率進行壓縮觀測，以保證重構(gòu)出的圖像質(zhì)量良好。

2.2 非參考波段的動態(tài)測量率觀測

2.2.1 提取邊緣區(qū)域

(3)

式中，x、y分別對應(yīng)圖像的橫軸方向和縱軸方向。接下來要計算平滑后的圖像I中每個像素點的梯度值和對應(yīng)的方向，先用一階微分算子分別得到點(i,j)兩個方向的偏導數(shù)Gx(i,j)和Gy(i,j)，那么圖像I中像素點的梯度幅值G(i,j)和梯度方向θ(i,j)分別為:

(4)

(5)

計算出全局的梯度仍然不能確定圖像邊緣，為達到細化梯度幅值圖像的目的，要進行極大值抑制處理，只保留梯度值的局部最大點。在圖像中以點(i,j)為中心的鄰域中沿梯度方向θ(i,j)進行插值，若點(i,j)處的梯度值大于相鄰插值點，則(i,j)為圖像梯度中的局部最大點，反之則為非局部極大值點，設(shè)置非局部極大值的點為0，最終得到細化后的圖像。

細化后的圖像是候選的圖像邊緣及邊緣線段，需要繼續(xù)篩選并連接出最終的圖像邊緣。Canny算法使用雙閾值法實現(xiàn)邊緣的檢測與連接。分別設(shè)置好一個高閾值和一個低閾值對細化后的圖像中記為局部最大值的點進行檢測，高于高閾值的點則為邊緣點，低于低閾值的點記為非邊緣點，處于兩閾值之間的點記為候選邊緣點，結(jié)合已確定為邊緣點的信息及連通性，若候選邊緣點臨近像素中已有邊緣點，則該候選邊緣點記為邊緣點，否則從候選邊緣點中去除。

用Canny算法提取出圖像邊緣后，需要劃分出圖像邊緣區(qū)域，確定哪些區(qū)域包含邊緣信息。為了更快速地劃分出圖像邊緣區(qū)域，本文使用圖像膨脹算法以將提取出的細化圖像邊緣增強擴展并劃分出圖像邊緣區(qū)域。具體操作為：將邊緣圖像與結(jié)構(gòu)元素進行卷積，本文使用大小為b×b的元素都為1的矩陣作為結(jié)構(gòu)元素，計算結(jié)構(gòu)元素在邊緣圖像中覆蓋區(qū)域像素點的最大值，并把這個最大值賦給結(jié)構(gòu)元素所覆蓋的區(qū)域，這樣就使邊緣膨脹加粗，最后得到記錄邊緣區(qū)域的圖像IE。

2.2.2 動態(tài)測量率設(shè)定方法

記錄邊緣區(qū)域的圖像IE是二值化圖像，值為1的像素點代表屬于邊緣區(qū)域，值為0的像素點代表屬于背景區(qū)域。在設(shè)定動態(tài)測量率之前，首先要對記錄邊緣區(qū)域的圖像分塊，文中利用對高光譜圖像分塊的方法對IE進行分塊，得到n個大小為B×B的不重疊塊，這樣分塊后可以使高光譜圖像中的圖像塊和IE中的塊互相對應(yīng)。接著計算每個圖像塊的像素點灰度值的總和S，第i塊的值記為Si,i=1,2,…,n，從而可以由S值知道每個圖像塊包含多少邊緣信息，S值高的圖像塊包含的邊緣信息多，S值低的圖像塊包含的邊緣信息少，S=0的塊不含邊緣信息即屬于背景區(qū)域。由于IE和高光譜分塊圖像的各個區(qū)域一一對應(yīng)，那么可知圖像塊xi對應(yīng)的邊緣信息值為Si。

在總測量值數(shù)目固定的條件下，需要對邊緣信息多的區(qū)域多采樣、背景區(qū)域少采樣從而達到提高圖像重構(gòu)質(zhì)量、優(yōu)化圖像細節(jié)的目的，即圖像塊的S值高，則該圖像塊的采樣率高，測量值的長度M大。本文建立一個線性模型將第i圖像塊的S值映射到測量值的長度，記Mi為第i個圖像塊要采集的測量值個數(shù)，那么映射關(guān)系如下:

Mi=kSi+m,

(6)

式中，k為待計算的比例系數(shù)，m為背景區(qū)域固定測量數(shù)目，m值的大小與總測量率、圖像塊像素數(shù)有關(guān)，設(shè)定要求為在不影響圖像整體重構(gòu)質(zhì)量的前提下使m值盡可能小。這個值可以在采集數(shù)據(jù)前根據(jù)經(jīng)驗值提前設(shè)定好。由于圖像包含了n個大小相同的圖像塊，那么每塊的測量值個數(shù)如式(7)所示：

M1=kS1+m,

M2=kS2+m,

?

Mn=kSn+m.

(7)

將上式等號左右兩邊分別全部相加，得到

(8)

(9)

將等式(9)中比例系數(shù)k的值代入邊緣信息與圖像塊測量值數(shù)目的映射關(guān)系式(6)中，則就可以利用邊緣信息得到每個圖像塊不同的自適應(yīng)測量值數(shù)目，則第i個圖像塊最終的動態(tài)測量值數(shù)目為Mi：

(10)

式中，round(·)表示對括號內(nèi)的計算結(jié)果四舍五入取整數(shù)。

3 利用邊緣信息進行動態(tài)測量流程

3.1 實施流程

對非參考波段的動態(tài)測量率自適應(yīng)采樣流程如圖2所示。

高光譜圖像利用邊緣信息進行動態(tài)測量的具體實施流程如下：

步驟1：將高光譜圖像的每個波段都分成n個大小為B×B的不重疊圖像塊；

步驟2：對參考波段的圖像按照固定采樣率采集，其采樣率比非參考波段設(shè)置的稍高一些，這樣可以保證參考波段的重構(gòu)質(zhì)量，從而使非參考波段自適應(yīng)分配更合理準確；

步驟3：將參考波段獲得的觀測數(shù)據(jù)用壓縮感知算法重構(gòu)出來，獲得該波段的圖像；

步驟4：對參考波段重構(gòu)出的圖像進行邊緣提取，并用圖像膨脹方法增強邊緣，劃分出圖像邊緣區(qū)域；

步驟5：提取的邊緣區(qū)域圖像也按照步驟1相同方法分成n個大小為B×B的不重疊圖像塊，計算每塊內(nèi)所有像素的灰度和，該值表示圖像塊含有的邊緣信息量，記為Si,i=1,2,…,n；

步驟6：設(shè)定背景區(qū)域(不含邊緣信息)的圖像塊的固定測量值數(shù)目m；

圖2 非參考波段動態(tài)測量率的自適應(yīng)采樣流程圖 Fig.2 Adaptive sampling flow chart of dynamic measurement rate in non-reference band

步驟7：根據(jù)公式(10)計算出每個圖像塊最終分配的自適應(yīng)測量值數(shù)目。

在高光譜圖像重構(gòu)端，接收到各個波段圖像的觀測數(shù)據(jù)后，采用相同的壓縮感知算法進行分塊重構(gòu)。

3.2 動態(tài)測量實施的可能性討論

在高光譜壓縮感知成像系統(tǒng)中，目前都是采用固定測量率的方式對圖像進行采集，尤其在分塊壓縮感知中，由于每個圖像塊都是同時采集時，必然導致直接獲取到的測量值數(shù)目相同。因此，動態(tài)分配測量值的壓縮感知算法就很難在系統(tǒng)中實現(xiàn)。如果必須要在高光譜壓縮感知系統(tǒng)中應(yīng)用動態(tài)測量算法，那么在系統(tǒng)中設(shè)定的采集測量率則不能小于所有圖像塊中最大測量率的值，即每個圖像塊采集次數(shù)為M≥max{Mi}，其中Mi是由式(10)計算出的每個圖像塊分配的測量值數(shù)目。這樣才能保證每個圖像塊都獲取到所需的動態(tài)測量值數(shù)量。

由于高光譜數(shù)據(jù)為三維數(shù)據(jù)，既有二維空間信息，又有很高的光譜分辨率，使得高光譜圖像數(shù)據(jù)量龐大。為了減小數(shù)據(jù)傳輸?shù)膲毫?，同時又保證圖像質(zhì)量，可以在采集一定數(shù)量的觀測數(shù)據(jù)后丟掉一些冗余信息來緩解數(shù)據(jù)傳輸壓力。在這種情況下，可以采用自適應(yīng)采樣率分配的算法來選擇保留不同圖像塊的觀測值數(shù)量。但由于傳統(tǒng)的自適應(yīng)測量率方法是根據(jù)圖像塊的方差大小分配測量率的，不同圖像紋理信息不同，極有可能造成圖像方差相差懸殊，導致不同圖像塊分配的測量值數(shù)目差異也非常大，最大分配數(shù)量甚至有可能超過圖像塊本身的像素數(shù)。而本文提出的利用邊緣信息動態(tài)測量高光譜圖像的方法中，對背景圖像塊也設(shè)置了固定測量值數(shù)目，并且按照圖像塊含有邊緣信息的多少來進行分配，從而可以避免某些圖像塊分配測量值數(shù)目過大的情況。

4 實驗結(jié)果與分析

本文對空間分辨率為256 pixel×256 pixel個像素，24個波段，波長為452～667 nm的Lego高光譜圖像進行仿真實驗。為了更準確地分配非參考波段的測量值數(shù)目，要選擇一個邊緣信息比較明顯的波段作為參考波段。這里選擇Lego高光譜圖像的第9波段作為參考波段優(yōu)先重構(gòu)。圖像分塊過程中，每塊的大小設(shè)置為8×8，參考波段的采樣率設(shè)置為0.5，用離散余弦基(DCT)作為稀疏矩陣，采用SpaRSA壓縮感知[8]算法獨立地重構(gòu)出各圖像塊。所有仿真實驗都在相同的計算平臺(Intel i5-4590 CPU 3.30 GHz /6.00 GB內(nèi)存)下進行。

圖3為參考波段采用固定測量率R=0.5時的重構(gòu)結(jié)果，使用的觀測矩陣為隨機01分布的矩陣，矩陣中0和1的比例為1∶1，SpaRSA算法中的參數(shù)值設(shè)為0.01。重構(gòu)出的圖像與原圖相比峰值信噪比(Peak Signal to Noise Radio,PSNR)為29.75 dB。

圖3 參考波段重構(gòu)圖像 Fig.3 Reconstructed image of key band

將重構(gòu)出的圖像作為先驗信息提取邊緣區(qū)域。下圖為對重構(gòu)后的參考波段圖像用Canny算子提取邊緣，并用5×5的結(jié)構(gòu)元素進行圖像邊緣膨脹后劃分出的圖像邊緣區(qū)域。

圖4 參考波段的邊緣區(qū)域劃分 Fig.4 Edge region extraction of key band

根據(jù)參考波段提取的邊緣信息，接下來要對非參考波段分配動態(tài)測量值數(shù)目。其中，有兩個參數(shù)可以根據(jù)經(jīng)驗值設(shè)置，分別為邊緣膨脹的寬度b，以及背景區(qū)域圖像塊(不含邊緣信息)的固定測量值m。下面以重構(gòu)第10波段為例，分析這兩個參數(shù)對動態(tài)測量算法重構(gòu)效果的影響。這里圖像總采樣率設(shè)為0.3，使用的觀測矩陣為隨機01分布的矩陣，矩陣中0和1的比例為1∶1，SpaRSA算法中的參數(shù)τ值設(shè)為0.01。由于觀測矩陣具有隨機性，每組實驗重構(gòu)圖像的PSNR值與重構(gòu)時間(s)均取5次測試的平均值。

總采樣率R=0.3，背景區(qū)域圖像塊m=5時，隨著邊緣膨脹寬度b的變化，第10波段的重構(gòu)質(zhì)量和重構(gòu)時間數(shù)值如表1所示。

表1 m=5時重構(gòu)結(jié)果

從表1可以看出，邊緣膨脹寬度b對重構(gòu)效果有影響，但是影響不很大，然而邊緣膨脹寬度b不能設(shè)置過大，超過5之后，PSNR有明顯下降趨勢，另外隨著b值增大，重構(gòu)時間逐漸增長。綜合表中結(jié)果來看，b=3是一個合適的邊緣膨脹寬度值。

接下來分析背景區(qū)域圖像塊的固定測量值m對重構(gòu)結(jié)果的影響。表2是總采樣率R=0.3，b=3時，第10波段的重構(gòu)質(zhì)量和重構(gòu)時間隨m值變化情況?？梢钥闯鲭S著固定測量值m的增大，圖像總體重構(gòu)效果稍稍變差，重構(gòu)時間增大。但是需要強調(diào)的一點是，在m=2時，許多圖像塊重構(gòu)失敗，所以固定測量值m值設(shè)定不能太小。m值設(shè)定過小時，雖然總體重構(gòu)質(zhì)量提高了，但是圖像會出現(xiàn)局部失真情況。另外，背景區(qū)域測量值m也與圖像總采樣數(shù)量有關(guān)，總采樣值高時m值也需要設(shè)置的相對高一些，避免圖像局部重構(gòu)效果過于失真。本文根據(jù)m=round(0.25·R·B2)來確定m的值,若計算出m<3時，m值設(shè)置為3。

基于上述參數(shù)設(shè)置分析結(jié)果，本文對Lego高

表2 b=3時重構(gòu)結(jié)果

光譜圖像中除第9波段外進行非參考波段重構(gòu)的仿真實驗。非參考波段的采樣率分別設(shè)置為0.1、0.2、0.3、0.4、0.5，邊緣膨脹寬度b=3，背景區(qū)域固定測量值m=round(0.25·R·B2)，表3與圖5為本文方法與傳統(tǒng)固定測量高光譜圖像的重構(gòu)結(jié)果比較，圖中結(jié)果為所有波段結(jié)果的平均值，并且每組實驗重復5次取平均。

從表3和圖5可以看出，本文提出的方法在總采樣率較小(R=0.1)時，PSNR提高了1 dB，但是重構(gòu)時間超過傳統(tǒng)固定測量值的方法。這是由于總采樣率過小的情況下，部分圖像塊無法獲得充足的測量次數(shù)，導致圖像質(zhì)量提升不顯著。但是在采樣率R大于0.1之后，利用邊緣信息動態(tài)測量的方法能顯著提升圖像質(zhì)量，PSNR增益為3～4 dB時，因為圖像邊緣區(qū)域獲得了足夠的測量以后，在圖像細節(jié)處圖像質(zhì)量改善明顯，隨著采樣率增大，本文提出方法的重構(gòu)時間也逐漸減小。

表3 非參考波段重構(gòu)PSNR(dB)對比

圖5 非參考波段重構(gòu)結(jié)果比較 Fig.5 Reconstruction results comparison for non-key bands

圖6給出了R=0.3，邊緣膨脹寬度b=3，背景區(qū)域固定測量值m=round(0.25·R·B2)時的測量值數(shù)目分配情況。其中單個圖像塊的測量值數(shù)目分配最大值為59，小于圖像塊像素數(shù)B2=64，所以在真實高光譜壓縮感知系統(tǒng)中可以實現(xiàn)。

圖6 測量值數(shù)目分配情況 Fig.6 Allocation of the number of measurements

圖7 第10波段重構(gòu)結(jié)果比較 Fig.7 Reconstruction results comparison for band 10

圖7是第10波段R=0.3時重構(gòu)效果對比圖，可以看到本文方法中玩偶的眼睛，衣服上的圖案等圖像細節(jié)處紋理清晰?？梢姳疚睦眠吘壭畔矸峙洳蓸勇实姆椒ㄊ怯行У?，提高了圖像的重構(gòu)質(zhì)量。

5 結(jié) 論

本文根據(jù)高光譜圖像譜間相關(guān)性很強的特點，提出了利用邊緣信息動態(tài)測量高光譜圖像的壓縮感知方法，通過根據(jù)圖像塊含有邊緣信息的多少自適應(yīng)分配采樣率，討論了動態(tài)測量方法在高光譜成像系統(tǒng)中實現(xiàn)的可行性。仿真結(jié)果表明，在總測量數(shù)目固定的條件下，采用本文提出的動態(tài)測量圖像塊方法重構(gòu)出的圖像比現(xiàn)有固定圖像塊采樣率采樣方法的質(zhì)量更高，PSNR比傳統(tǒng)方法高出1～4 dB圖像細節(jié)處更清晰。下一步的工作重點將研究如何進一步合理分配測量值數(shù)目并將其成功應(yīng)用于真實高光譜壓縮感知系統(tǒng)。