陸基高光譜成像用于核驅(qū)動模型測量與地物分類

王強輝， 沈?qū)W舉， 周 冰， 華文深， 應家駒， 趙佳樂

陸軍工程大學石家莊校區(qū)電子與光學工程系， 河北 石家莊 050003

引 言

高光譜成像是一種先進的圖像獲取技術(shù)， 利用成像光譜儀獲取地物的空間和光譜信息， 光譜分辨率高， 波段寬度通常在10 nm以下， 利用這一特點可以區(qū)分出不同的地物類型[1]。 高光譜地物分類是指利用高光譜圖像中的空間和光譜信息， 對圖像中的未知像元或未知地物進行類別標注， 判斷所屬類別[2]。 高光譜圖像中不同的地物類型通常具有不同的幾何空間信息和光譜信息， 根據(jù)分類前是否已獲取地物的光譜信息， 地物分類可分為監(jiān)督分類和非監(jiān)督分類。

傳統(tǒng)的高光譜成像方式多采用空基或者天基的方式， 這兩種成像方式成像高度高， 偵察時間固定， 探測方向基本垂直于地面， 拍攝對象一般為大尺度的目標。 然而， 隨著成像光譜儀的發(fā)展， 特別是光譜儀的不斷小型化， 陸基應用成為可能。 陸基成像是指使用手持或者小型無人機為載體的成像方式， 其成像時間任意， 探測角度任意， 與傳統(tǒng)的高光譜成像方式有著較大的差別。 陸基成像條件下的地物光譜受入射天頂角、 光源與探測器相對方位角、 探測天頂角等因素影響顯著， 同種地物在不同成像條件下光譜數(shù)據(jù)往往發(fā)生一定的變化， 這種變化會對后續(xù)圖像的處理比如地物分類產(chǎn)生影響[3]。

二向反射分布函數(shù)(bidirectional reflectance distribution function， BRDF)模型[4]建立了方向反射與地物參數(shù)的關(guān)系， 反映了地物對太陽輻射的反射能力。 目前比較成熟的BRDF模型有物理模型、 經(jīng)驗模型和半經(jīng)驗模型等。 其中半經(jīng)驗模型具有簡單、 可操作性強的優(yōu)點， 并且外推可以求得任意成像條件下的地物二向反射特征， 反演速度快， 適用于陸基條件下高光譜成像的研究[5]。 核驅(qū)動模型是半經(jīng)驗模型的一種， 具有擬合能力強的優(yōu)點。 李小文等[6]證明了半經(jīng)驗核驅(qū)動模型的擬合能力與反演能力。 程晨等[7]研究了多種植物在不同波段范圍的擬合效果， 結(jié)果表明半經(jīng)驗核驅(qū)動模型具有很好的反演能力。

利用三組不同成像條件的地物數(shù)據(jù)測量了各地物的散射系數(shù)， 并根據(jù)散射系數(shù)對第四組數(shù)據(jù)進行擬合證實了核驅(qū)動模型的準確性。 通過對散射系數(shù)的表達式f(λ)可以看出散射系數(shù)只與地物類型和波長有關(guān)， 與成像條件無關(guān)， 實驗結(jié)果發(fā)現(xiàn)不同地物的散射系數(shù)有著較大差異， 因此可以利用這種差異對不同地物進行分類， 并用投影、 距離、 信息量三種相似性度量指標證明了分類的有效性。 這種分類方法不受成像條件的影響， 在獲得地物的散射系數(shù)向量后可以對任意成像條件下的地物分類， 具有重要的應用價值。

1 BRDF模型及陸基測量原理

自然地物絕大部分是非朗伯體， 一般使用BRDF模型描述其各向異性。 BRDF模型定義為反射幅亮度與入射輻照度之比， 反映了地物對入射光線的反射情況。

1.1 BRDF模型

BRDF定義為光源在入射點附近面元上的反射幅亮度與入射輻照度之比， 其計算表達式為

(1)

式(1)中，φi,θi,φr,θr分別表示入射方位角、 入射天頂角、 反射方位角和反射天頂角。 dE(φi,θi)表示光源在入射點附近面元上的入射輻照度， dL(φi,θi;φr,θr)為相應的入射點附近面元的反射輻亮度。

1.2 半經(jīng)驗核驅(qū)動模型

根據(jù)計算角度的不同， BRDF模型可分為物理模型、 經(jīng)驗模型和半經(jīng)驗模型。 半經(jīng)驗模型具有快速簡潔、 適應性強的特點， 應用更為廣泛。 核驅(qū)動模型是半經(jīng)驗模型的一種， 其中Rossthick-LiSpareR模型和Rossthick-LiTransitN模型最具代表性。 本文選擇的半經(jīng)驗核驅(qū)動模型為Rossthick-LiSpareR模型， 其線性組合各向同性散射核、 體散射核、 幾何光學散射核來描述地物的二向反射特性， 表達式如式(2)所示

Rs(θi,θr,σ,λ)=fiso(λ)+fvol(λ)Kvol(θi,θr,σ)+

fgeo(λ)Kgeo(θi,θr,σ)

(2)

該模型將二向反射特性分解為三個部分， 分別為各向同性散射、 體散射以及幾何光學散射， 其中Rs(θi,θr,σ,λ)表示二向反射率，θi表示入射天頂角，θr表示觀測天頂角，σ表示相對方位角，Kvol表示體散射核，Kgeo表示幾何光學核， 其體散射核是Rossthick(羅斯厚層)核， 幾何光學核是LiSpareR(李氏稀疏互易)核。fiso，fvol和fgeo為三個只與波長和地物類型有關(guān)的散射系數(shù)， 分別表示均勻散射、 體散射和幾何光學散射所占比例。

RossThick核的表達式如式(3)所示

(3)

式(3)中，δ是相位角， cosδ=cosθicosθr+sinθisinθrsinσ。

LiSpareR核的表達式為

Kgeo(θi,θr,σ)=Ο(θi,θr,σ)-secθi-

(4)

式(4)中

(5)

(6)

(7)

體散射核和幾何光學核可以通過探測天頂角、 相對方位角以及入射天頂角進行計算。 在計算二向反射率的時候， 由于兩個核的積分與常系數(shù)無關(guān)， 可以先計算出兩個核的值， 利用至少三組不同條件下的地物二向反射率與核的值可求得三個散射系數(shù)的值， 進而利用三個散射系數(shù)進行外推可以求得任意入射天頂角、 相對方位角以及觀測天頂角條件下的二向反射率。 入射天頂角、 相對方位角以及觀測天頂角分布如圖1所示。

圖1 入射天頂角、 相對方位角以及觀測天頂角分布Fig.1 Schemadic of Incident zenith angle, relative azimuth angle, and observation zenith angle

1.3 陸基測量原理

實驗在陸基成像的條件下進行， 陸基測量存在天空漫反射光， 對實驗會產(chǎn)生影響， 為了更加準確地測得地物的二向反射特性， 需要同時考慮到來自太陽的直射光和天空中的漫反射光。 測量得到的反射率為兩者之和， 地物二向反射率為測量得到的總的反射率與測得的漫反射反射率之差， 見式(8)

R(θi,φi;θr,φr)=K1Rs(θi,φi;θr,φr)+K2RD(θr,φr)

(8)

式(8)中

(9)

(10)

式中，θi為太陽的天頂角，φi為太陽的方位角，θr為觀察儀器的天頂角，φr為觀察儀器的方位角，ID為天空漫反射光照射地物的輻射強度，Is(θi,φi)為太陽直射光照射地物的輻射強度，I(θi,φi)為太陽直射光和漫反射光的總輻射強度，RD(θr,φr)為漫入射的半球一定方向的反射率，Rs(θi,φi;θr,φr)為太陽直射光照射下的二向反射率，R(θi,φi;θr,φr)為陸基條件下測量出的總反射率。

由于實地情況較為復雜， 我們認為目標與標準板的測量值之比就是反射率之比。 目標反射率的計算式為

(11)

式(11)中，ρs(λ)為標準板的反射率，V(λ)為目標物的儀器測量值，Vs(λ)為標準板的儀器測量值。 通過式(8)即可得到只有太陽直射條件下的地物二向反射率Rs， 即

(12)

2 實驗部分

2.1 儀器及參數(shù)

選用的儀器是基于聲光可調(diào)諧濾波器(Acousto-Optic Tunable Filer, AOTF)的HSI-300型成像光譜儀， 光譜分辨率為4 nm， 空間分辨率為1 002×1 002像元， 在449～801 nm的波段范圍內(nèi)獲取89幅不同波段的灰度圖像， 每一幅圖像都記錄了在不同波段下地物的輻射強度值。 實驗地點選擇在河北省石家莊市某地， 地理坐標為北緯38°27′， 東經(jīng)114°30′， 時間為2021年4月27日， 天氣晴朗， 實驗時氣溫為20～24 ℃， 分別在11:10， 11:40， 12:10和14:10進行四組實驗。 具體實驗環(huán)境如表1所示， 在不同條件下共得到8組數(shù)據(jù)。 數(shù)據(jù)1-1， 2-1， 3-1和4-1為在太陽光直射、 無遮擋條件下的情況， 數(shù)據(jù)1-2， 2-2， 3-2和4-2為用擋板擋住太陽光照射、 有遮擋條件下的情況, 如表1所示。

表1 實驗條件記錄Table 1 Experimental conditions record

選取的研究對象是常見迷彩及偽裝板， 迷彩主要包括迷彩雨衣B、 迷彩水壺C、 迷彩鞋D、 荒漠迷彩帽E、 叢林迷彩帽F以及假草皮G； 偽裝板主要有叢林綠(深)H、 叢林綠(淺)J、 荒漠色(深)I、 荒漠色(淺)K四種。 地物擺放示意圖如圖3所示， 圖4為各種地物在數(shù)據(jù)4-1條件下的光譜曲線。

圖2 有無遮擋條件下測量方法示意圖(a)： 無遮擋； (b)： 有遮擋Fig.2 Schematic diagram of the measurement with or without obstruction(a)： Without obstruction； (b): With obstruction

圖3 地物擺放示意圖Fig.3 Schematic diagram of the placement of ground features

圖4 地物光譜曲線Fig.4 Spectra of ground features

為了更加準確測定地物反射特性， 采用白板法對地物光譜數(shù)據(jù)進行輻射定標。 采用的標準白板A為經(jīng)過計量標定的聚四氟乙烯板(PTFE), PTFE材料各方向的反射特性比較均勻， 可以近似的看為一個朗伯體。 PTFE板反射率θ=0.98。 通過同時拍攝白板和地物可以求得地物的光譜反射率， 其計算公式為

(13)

2.2 方法

按照以下幾個步驟進行實驗， 將實驗人員分為地面人員和探測人員兩組， 地面人員主要負責對目標地物和白板進行設(shè)置， 并適時用擋板擋住太陽直射光， 配合探測人員進行實驗； 探測人員主要操作成像光譜儀對目標地物和白板進行成像， 并指揮地面人員進行操作。

2.2.1 測量K2和K1

(1) 在自然太陽直射條件下， 地面人員將白板放置在空曠草地上， 操作人員利用成像光譜儀從七樓窗臺對地面的白板進行成像， 測量得到不同光譜段的輻射強度I。

(2) 地面人員迅速用擋板遮住太陽光使陰影完全蓋住白板， 再利用成像光譜儀對白板進行測量得到ID， 此時處于無太陽直射光， 僅存在天空中的漫反射光， 由于遮光過程較短， 此時可將入射天頂角、 探測天頂角、 相對方位角等信息視為不變；

(3) 利用matlab對公式K2=ID/I(θi,φi)求解計算得到K2， 再由K1=1-K2得到K1。

2.2.2 測量R和RD并求解RS

(1) 在11:10， 11:40， 12:10和14:10四個時間點， 利用成像光譜儀迅速測量目標地物和白板的輻射值， HSI-300型成像光譜儀視場較大， 可同時對目標地物和標準板進行成像。 入射天頂角、 探測天頂角、 相對方位角在每次測量都認為保持不變， 通過目標地物的輻射值與標準漫反射板的輻射值相比得到R(θi,φi;θr,φr)， 地面人員記錄每組實驗的入射天頂角θi、 探測天頂角θr以及太陽光入射方向與探測方向的相對方位角σ， 入射天頂角θi由石家莊地區(qū)經(jīng)緯度、 海拔以及實驗的時間查詢得到， 例如2021年4月27日11:10進行實驗時太陽高度角約為27°； 利用電子羅盤和經(jīng)緯儀測量得到探測天頂角θr， 通過投影法測量得到相對方位角σ；

(2) 地面人員用擋板擋住太陽直射光， 在只存在天空漫反射光時， 探測人員再次測量目標地物和標準板的輻射值， 利用matlab計算目標地物和標準板的輻射值的比值得到該方向的漫入射的半球一定方向反射率RD(θr,φr)；

(3) 根據(jù)式(2)—式(8)計算得到Rs(θi,φi;θr,φr)， 即只有太陽直射光時地物的反射率。

2.2.3 計算散射系數(shù)

(1) 利用測量得到的入射天頂角θi、 探測天頂角θr、 相對方位角σ求解Kvol，Kgeo；

(2) 將前三組實驗測量計算得到的Kvol，Kgeo，Rs(θi,φi;θr,φr)代入模型公式， 利用線性擬合的方法分別得到地物的fiso，fvol，fgeo值。 利用matlab對實驗數(shù)據(jù)進行處理， 得到迷彩雨衣、 迷彩水壺、 迷彩鞋、 荒漠迷彩帽、 叢林迷彩帽、假草皮和叢林綠(深)、 叢林綠(淺)、 荒漠色(深)、 荒漠色(淺)四種偽裝板以及草地的散射系數(shù)分布。

2.3 數(shù)據(jù)處理

為了驗證Rossthick-LiSpareR模型的擬合能力， 利用前三組數(shù)據(jù)求解散射系數(shù)， 將第四組的觀測條件， 即已測的入射天頂角， 探測天頂角、 相對方位角代入求解得到Kvol和Kgeo， 再將Kvol和Kgeo與散射系數(shù)代入得到擬合的光譜RS， 利用第四組數(shù)據(jù)的實測值對擬合的光譜進行檢驗， 如圖6所示。

通過圖5和圖6可以發(fā)現(xiàn)， 散射系數(shù)是反映地物類別的特征， 其中fiso和fvol比較相似， 并且散射系數(shù)中出現(xiàn)了負值情況， 這是因為地物的均勻散射、 體散射、 幾何光學散射并不是完全獨立的， 與地物的類型、 空間分布等特點密切相關(guān)。 并且本實驗證明了半經(jīng)驗核驅(qū)動模型的良好擬合能力， 對于絕大多數(shù)地物擬合光譜與實測光譜有著極高的相似度。 但也存在一些波段擬合效果稍差， 比如假草皮與迷彩鞋在650～800 nm光譜差異相對較大， 可能是由于BRDF模型本身對不同分布地物的擬合性能差異或者儀器本身偏差、 天氣等原因造成。

圖5 不同地物散射權(quán)系數(shù)分布Fig.5 Scattering coefficients of different ground objects

圖6 不同地物擬合光譜與實測光譜Fig.6 Fitted and measured spectra of different ground features

3 結(jié)果與討論

通過圖5可以發(fā)現(xiàn)， 不同類型地物的散射系數(shù)有較大差別， 散射系數(shù)的分布反映了地物的類型， 并且與成像條件或者探測角度等條件無關(guān)， 不受地物二向反射特征的影響， 因此可以根據(jù)地物之間散射系數(shù)的相似性對未知地物進行分類。 分類的思路是首先對未知某一地物求解fiso，fvol和fgeo值， 將其fiso，fvol和fgeo值與已知地物的fiso，fvol和fgeo值進行對比， 比較二者之間的相似性， 相似度越高則越有可能是該已知地物類別， 從而對未知地物進行分類。 選擇兩組數(shù)據(jù)對分類效果進行驗證， 數(shù)據(jù)1為截取的大小為200×200像元空間范圍內(nèi)的一未知地物， 如紅框區(qū)域所示， 選取示意圖與其散射系數(shù)如圖7所示。

圖7 未知地物選取示意圖與其散射系數(shù)Fig.7 Schematic diagram of the selection of unknown features and their scattering coefficients

數(shù)據(jù)2為一種合成地物， 迷彩雨衣、 迷彩水壺、 迷彩鞋、 荒漠迷彩帽、 叢林迷彩帽、 假草皮和叢林綠(深)、 叢林綠(淺)、 荒漠色(深)、 荒漠色(淺)四種偽裝板、 草地各按3%， 3%， 3%， 3%， 3%， 3%， 3%， 3%， 3%， 3%和70%的比例合成， 合成地物的光譜曲線及散射系數(shù)如圖8所示。

圖8 合成地物光譜與其散射權(quán)系數(shù)Fig.8 Spectra of synthetic ground objects and their scattering coefficients

從兩組散射系數(shù)數(shù)據(jù)可以定性地看出， 其散射系數(shù)與草地較為相似， 因此該地物大概率為草地。 為了定量地衡量其與各地物的相似性， 分別用投影、 距離、 信息量三種相似性度量指標來衡量曲線間相似性。

3.1 基于投影的相似性度量方法

光譜角制圖(spectral angle metric， SAM)法[8]是一種經(jīng)典的相似性度量方法， 反映了光譜曲線間形狀的差異， 其值對應的是兩曲線的余弦夾角。 不同類型地物間散射系數(shù)的差異同樣可以采用光譜角制圖法進行分析。 設(shè)某地物和其他已知地物的某個fiso，fvol和fgeo向量分別為M和Mc， 定義光譜角距離為

(14)

式(14)中，MT為M的轉(zhuǎn)置。 SA(M,Mc)是向量M,Mc間的廣義夾角， 其值越小說明其相似性越好。 為了分析兩地物散射系數(shù)的綜合相似性， 將fiso，fvol和fgeo值間的光譜角距離以一定的權(quán)值相加得到綜合光譜角相似度dSAM

dSAM=ω1×SA(M,Mc)fiso+ω2×SA(M,Mc)fgeo+

ω3×SA(M,Mc)fvol

(15)

式(15)中， (ω1,ω2,ω3)為權(quán)值向量， 各取1/3。dSAM為綜合光譜角相似度，dSAM越小說明相似度越大， 極有可能屬于同種地物， 相反dSAM越大說明相似度越小， 屬于同種地物概率越低。

3.2 基于距離的相似性度量方法

均方根誤差(root mean square error， RMSE)法[9]是一種經(jīng)典的相似性度量方法， 反映了兩地物光譜矢量大小的差異， 因此均方根誤差可以運用于地物的相似度評價。 RMSE的公式如式(16)

(16)

式(16)中，n為光譜波段數(shù)，ED為歐氏距離， 根號中的系數(shù)1/n去除了光譜矢量大小對光譜維數(shù)的相關(guān)性， 從而RMSE表示兩個光譜矢量間的平均距離， 值越小表示其矢量距離越小。 為分析兩地物間散射系數(shù)的綜合相似性， 將fiso，fvol和fgeo三值的矢量距離以一定的權(quán)值相加得到綜合矢量距離相似度dRMSE

dRMSE=ω1×RMSE(M,Mc)fiso+ω2×RMSE(M,Mc)fgeo+

ω3×RMSE(M,Mc)fvol

(17)

式(17)中， (ω1,ω2,ω3)為權(quán)值向量， 在此各取1/3。dRMSE越小說明相似性越大， 兩地物屬于同種地物可能性越大， 相反dRMSE越大說明相似性越小， 屬于同種地物可能性越小。

3.3 基于信息量的相似性度量方法

熵描述了系統(tǒng)的混亂程度， 信息熵[10]描述了信號的不確定程度， 反映了信號信息量的大小。 其公式定義為

(18)

互信息[11]用于描述兩個系統(tǒng)間的相關(guān)性， 即一個系統(tǒng)中包含另一個系統(tǒng)中信息量的大小。 兩個信號的互信息(mutual information)定義為

(19)

將對數(shù)項進一步分解

-logp(x)-logp(y)-(-logp(x,y))

(20)

式(20)表明， 兩信號的互信息為兩信號的信息量之和減去兩信號同時發(fā)生的信息量后的剩余信息量， 反映了兩信號單獨發(fā)生重復信息量的大小。 互信息越大， 重復信息量越大， 兩信號相關(guān)性越大， 信號越相似。 將地物之間的互信息運用于相似度評價。 為分析兩地物間散射系數(shù)的綜合相似性， 將fiso，fvol和fgeo三者值的互信息以一定的權(quán)值相加得到綜合信息量相似度dMI

dMI=ω1×MI(M,Mc)fiso+ω2×MI(M,Mc)fgeo+

ω3×MI(M,Mc)fvol

(21)

式(21)中， (ω1,ω2,ω3)為權(quán)值向量， 在此各取1/3。dMI越大則相似性越大， 兩地物屬于同種地物可能性越大， 反之則相似性越小， 屬于同種地物可能性越小。

3.4 相似性計算與分類

分別計算兩組數(shù)據(jù)與各已知地物的綜合光譜角相似度dSAM、 綜合矢量距離相似度dRMSE與綜合信息量相似度dMI， 結(jié)果如表2所示。

表2 相似度計算Table 2 Similarity calculation

通過分析發(fā)現(xiàn)， 兩組數(shù)據(jù)中的地物與L(草地)的綜合光譜角相似度最小、 綜合矢量距離相似度最小并且綜合信息量相似度最大， 所以該地物屬于草地的可能性最大。 在選取數(shù)據(jù)1未知地物時， 該處地物分布確實屬于草地， 因此該組實驗證明了利用散射系數(shù)進行地物分類的有效性。 數(shù)據(jù)2是一組主要由草地構(gòu)成的合成數(shù)據(jù)， 實驗證明了在混合了多種不同類型地物后仍然可以根據(jù)散射系數(shù)將合成數(shù)據(jù)中的主要地物類型準確分類， 具有重大的實踐意義。

4 結(jié) 論

通過實驗的方法對多種地物的BRDF模型各參數(shù)進行了測量與求解， 證實了半經(jīng)驗核驅(qū)動模型良好的擬合效果。 陸基條件下高光譜成像地物光譜受成像條件和二向反射特性影響顯著， 光譜不固定， 在未知具體實驗條件或探測條件下無法直接利用光譜曲線對地物類別進行判別。 針對這一問題， 考慮到散射系數(shù)是地物的固有屬性， 與成像條件或者探測條件無關(guān)， 提出了一種利用散射系數(shù)進行分類的方法， 并通過綜合光譜角相似度、 綜合矢量距離相似度和綜合信息量相似度對該方法進行了驗證， 結(jié)果表明， 該方法克服了陸基條件下地物二向反射特性導致光譜不固定的缺點， 將這一過程進行外推， 對于其他任意未知地物， 計算其散射系數(shù)并與已知地物的散射系數(shù)進行對比， 求解相似度就可以對其進行分類， 具有重要的應用價值。