李 力, 唐曉霏, 劉 雯, 陳 城

(1. 湖北省氣象信息與技術(shù)保障中心, 武漢 430074； 2. 湖北省測(cè)繪工程院, 武漢 430074)

無(wú)人機(jī)遙感技術(shù)目前已廣泛應(yīng)用于氣象環(huán)境監(jiān)測(cè)、 國(guó)土測(cè)繪勘察、 水文地質(zhì)救災(zāi)等領(lǐng)域, 但由于無(wú)人機(jī)主要為中低空飛行, 而水汽云霧和霧霾天氣都會(huì)影響可見(jiàn)光和紅外線在大氣中的傳播, 因此使其獲取的圖像模糊有霧[1]. Oakley等[2-3]首次提出了有霧圖像的退化處理, 目前對(duì)單幅圖像的去霧算法研究已有很多結(jié)果[4-5]. He等[6]采用暗通道先驗(yàn)(dark channel prior, DCP)原理量化大氣透射率, 獲取大氣透射率圖層(transmission map, TM), 再用軟摳圖方法得到相對(duì)準(zhǔn)確的大氣透射率, 建立物理模型恢復(fù)清晰影像. 但該方法對(duì)場(chǎng)景中存在大量白色區(qū)域的情況去霧效果不佳, 相關(guān)算法復(fù)雜程度較高[7-8]、 計(jì)算速度較慢, 難以適應(yīng)航拍圖像數(shù)據(jù)量大、 及時(shí)性強(qiáng)等特點(diǎn).

通過(guò)分析基于暗通道先驗(yàn)的圖像去霧算法, 發(fā)現(xiàn)運(yùn)行中耗時(shí)較多的兩個(gè)環(huán)節(jié)是計(jì)算暗通道圖像和使用引導(dǎo)濾波精細(xì)化透射率圖層. 本文針對(duì)引導(dǎo)濾波精細(xì)化大氣透射率圖層時(shí)間較長(zhǎng)的問(wèn)題, 考慮采用與濾波窗口無(wú)關(guān)、 算法時(shí)間復(fù)雜度為O(1)的中值濾波精細(xì)化透射率圖層. 改進(jìn)算法通過(guò)計(jì)算機(jī)編碼, 利用無(wú)參考的客觀圖像質(zhì)量評(píng)價(jià)方法進(jìn)行量化評(píng)價(jià), 將信息熵、 邊緣強(qiáng)度和方差作為評(píng)價(jià)圖像質(zhì)量的3個(gè)指標(biāo), 以在去霧圖像質(zhì)量和后期處理應(yīng)用影響較小的條件下節(jié)省計(jì)算時(shí)間, 提高航拍圖像暗通道先驗(yàn)去霧算法的處理速度.

1 暗通道先驗(yàn)圖像去霧方法

1.1 有霧圖像退化模型

在計(jì)算機(jī)視覺(jué)圖形學(xué)領(lǐng)域, 有霧圖像的退化模型[9]表示為

I(x)=J(x)t(x)+A(1-t(x)),

(1)

其中:I(x)為傳感器接收到的場(chǎng)景信號(hào), 即輸入的有霧影像； 場(chǎng)景輻射J(x)為信號(hào)處理后的清晰圖像；t(x)為大氣透射率,J(x)t(x)稱(chēng)為直接衰減項(xiàng), 用于量化場(chǎng)景輻射和傳播中的信號(hào)損失；A為環(huán)境光照強(qiáng)度,A(1-t(x))稱(chēng)為環(huán)境光照, 表示圖像接收的大氣散射光強(qiáng)[10], 是引起色彩偏移和云霧效果的直接原因.

1.2 暗通道先驗(yàn)理論

基于統(tǒng)計(jì)大量清晰圖像得到的暗通道先驗(yàn)理論, 是指大部分不含天空優(yōu)質(zhì)圖像的所有像素在R,G,B 3個(gè)通道中最少存在一個(gè)顏色通道灰度值相當(dāng)?shù)蜕踔邻吔?, 即在一定的微小區(qū)域內(nèi)最小輻射強(qiáng)度值極低. 對(duì)于一副圖像J, 可定義

(2)

其中Ω(x)是中心為像素x的任一微小區(qū)域,y是該區(qū)域內(nèi)任一像素,Jc(y)為像素y在c通道的灰度值. 在不含天空區(qū)域的其他范圍內(nèi), 清晰圖像的Jdark通道灰度值均極小甚至趨近于0, 圖像J的暗通道值即為Jdark, 該統(tǒng)計(jì)理論稱(chēng)為暗通道先驗(yàn).

1.3 暗通道先驗(yàn)圖像去霧

去霧目標(biāo)是將傳感器接收的有霧圖像利用獲得的有用信息通過(guò)去霧還原出清晰圖像. 暗通道先驗(yàn)圖像去霧方法先根據(jù)暗通道原理獲取先驗(yàn)知識(shí), 再用有霧圖像退化模型實(shí)現(xiàn)圖像去霧.

1.3.1 估算環(huán)境光照強(qiáng)度

暗通道圖像中灰度值越高的區(qū)域云霧越厚, 在輸入圖像中位于該部分區(qū)域的像素灰度值越接近于環(huán)境光照強(qiáng)度. 估算環(huán)境光照強(qiáng)度首先需找到暗通道圖像中灰度值最高并占圖像總像素?cái)?shù)量0.1%的像素點(diǎn), 記錄其對(duì)應(yīng)的坐標(biāo)索引, 然后根據(jù)坐標(biāo)索引在輸入的有霧圖像中找到對(duì)應(yīng)像素點(diǎn), 計(jì)算有霧圖像中對(duì)應(yīng)像素點(diǎn)的灰度平均值作為環(huán)境光照強(qiáng)度A.

1.3.2 估算大氣透射率

(3)

按下式再次運(yùn)算3個(gè)顏色通道的最小值, 可得以像素x為中心的濾波窗口內(nèi)的灰度最小值:

(4)

根據(jù)暗通道先驗(yàn)原理和式(2)可得

(5)

于是, 由式(4)和式(5)可計(jì)算出大氣透射率為

(6)

1.3.3 圖像去霧處理

暗通道先驗(yàn)條件用于量化云霧厚度和全部像素的輻射還原量, 恢復(fù)出優(yōu)質(zhì)清晰的圖像. 通過(guò)式(1)的有霧圖像退化模型和環(huán)境光照強(qiáng)度A、 大氣透射率t(x), 可進(jìn)行單幅圖像去霧處理, 表示為

(7)

2 精細(xì)化大氣透射率圖層

通過(guò)對(duì)暗通道先驗(yàn)算法的剖析和代碼實(shí)驗(yàn)發(fā)現(xiàn), 計(jì)算暗通道圖像和使用引導(dǎo)濾波精細(xì)化透射率圖層時(shí)運(yùn)行耗時(shí)較多, 下面通過(guò)改進(jìn)大氣透射率圖層精細(xì)化處理方法提高去霧運(yùn)算的速度.

2.1 大氣透射率圖層精細(xì)化處理

根據(jù)上述暗通道圖像得到的大氣透射率仍不能滿(mǎn)足去霧精度要求, 此時(shí)暗通道圖像只能大致展現(xiàn)云霧的厚薄情況. 由于在求取暗通道過(guò)程中進(jìn)行了濾波窗口操作, 若用式(7)直接去霧會(huì)產(chǎn)生光暈現(xiàn)象. 在暗通道圖層中, 任意像素暗通道值都是以該像素為中心的n×n窗口區(qū)域內(nèi)所有像素在R,G,B通道的最低灰度值, 若圖層中有一個(gè)像素最低灰度值很小, 則以該像素為中心的n×n區(qū)域內(nèi)所有像素暗通道值都將很小, 導(dǎo)致區(qū)域灰度偏低, 形成擴(kuò)張現(xiàn)象, 因此去霧前需精細(xì)化處理大氣透射率圖層.

2.2 引導(dǎo)濾波精細(xì)化透射率圖層

引導(dǎo)濾波(guided filter, GF)[11-12]的處理效果與雙邊濾波相似, 但算法時(shí)間復(fù)雜度為O(n). 使用引導(dǎo)濾波的輸出圖像與引導(dǎo)圖像是局部線性關(guān)系, 處理圖像的平滑效果和邊緣保護(hù)效果更顯著. 引導(dǎo)濾波使輸出圖像與引導(dǎo)圖像具有相似的紋理結(jié)構(gòu)和場(chǎng)景信息[13-14], 在圖像去霧算法中可用于精細(xì)化大氣透射率圖層, 處理過(guò)程如圖1所示.

圖1 引導(dǎo)濾波示意圖Fig.1 Schematic diagram of guided filtering

線性平移變量濾波(linear translation-variant filtering, LTVF)中包括輸入圖像p、 輸出圖像q和引導(dǎo)圖像I. 輸入圖像與引導(dǎo)圖像是預(yù)先設(shè)定的兩幅相同或相異圖像, 特定像素i的濾波過(guò)程為

(8)

其中i,j表示像素坐標(biāo)索引,Wij為引導(dǎo)圖像相關(guān)的濾波核, 其不依賴(lài)于圖像p, 引導(dǎo)濾波器與輸入圖像p線性相關(guān)[15-16].

設(shè)引導(dǎo)濾波是由輸出圖像q與引導(dǎo)圖像I構(gòu)成的線性相關(guān)模型, 表示為

qi=akIi+bk, ?i∈wk,

(9)

其中線性系數(shù)(ak,bk)在以像素k為中心的濾波窗口區(qū)域wk內(nèi)是定值. 由式(9)可知, 輸出圖像q的梯度隨引導(dǎo)圖像I的梯度變化而改變, 引導(dǎo)濾波可將梯度和紋理信息傳遞給輸出圖像, 使兩種圖像的結(jié)構(gòu)和場(chǎng)景信息非常接近, 因此引導(dǎo)濾波可用于精細(xì)化大氣透射率圖層、 圖像自動(dòng)分割等領(lǐng)域. 為計(jì)算線性系數(shù)(ak,bk), 假設(shè)輸出圖像q是輸入圖像p減去噪聲部分, 表示為

qi=pi-ni,

(10)

式中ni表示噪聲.

在線性模型式(9)的基礎(chǔ)上, 要使輸入圖像p中噪聲占比最小, 需要輸入與輸出圖像間的差異最小. 用期望值E(ak,bk)量化以像素k為中心的濾波窗口區(qū)域內(nèi)輸入圖像p與輸出圖像q之間的差異, 表示為

(11)

其中ε為正則化參數(shù), 用于確定ak在式中所占比率. 如需濾波窗口區(qū)域wk內(nèi)噪聲占比最小, 則求解使式(11)取最小值的線性參數(shù)對(duì)(ak,bk)為

(12)

(13)

(14)

將線性系數(shù)(ak,bk)通過(guò)式(9)計(jì)算得到輸出圖像q, 再用式(8)可得出引導(dǎo)濾波內(nèi)核為

(15)

在圖像去霧過(guò)程中, 使用引導(dǎo)濾波對(duì)計(jì)算出的大氣透射率圖層精細(xì)化處理前后對(duì)比效果如圖2所示.

圖2 大氣透射率圖層經(jīng)引導(dǎo)濾波處理前后對(duì)比Fig.2 Comparison of atmospheric transmittance layers before and after guided filtering

2.3 中值濾波精細(xì)化透射率圖層

暗通道先驗(yàn)是在一個(gè)微小的局部區(qū)域內(nèi), 需要很低的最小輻射強(qiáng)度[17-18], 暗通道濾波用濾波方式找到其窗口區(qū)域內(nèi)輻射強(qiáng)度的最小值[19-21]. 將暗通道濾波視為二維黑白圖像的最小值濾波向三維空間彩色圖像的擴(kuò)展, 三維空間指R,G,B顏色空間. 暗通道濾波要找到濾波窗口區(qū)域內(nèi)3個(gè)顏色空間中輻射強(qiáng)度最小的像素, 并將該像素的輻射強(qiáng)度作為暗通道值. 為避免光暈現(xiàn)象, 需選取一個(gè)適合的濾波器, 在實(shí)現(xiàn)圖像平滑處理的同時(shí)保護(hù)其邊緣信息. 中值濾波器相比引導(dǎo)濾波更易實(shí)現(xiàn), 且算法時(shí)間復(fù)雜度較低, 能節(jié)約圖像去霧時(shí)間. 由于兩個(gè)直方圖累加是一個(gè)O(1)操作, 使得對(duì)每增加一個(gè)像素, 直方圖只需累加一個(gè)固定的數(shù), 從而獲得O(1)的復(fù)雜度. 實(shí)現(xiàn)算法時(shí)間復(fù)雜度為O(1)的最小值濾波過(guò)程如下.

1) 初始化核直方圖. 若圖像有n列, 每個(gè)直方圖對(duì)應(yīng)于圖像中的一列, 則需處理n個(gè)直方圖, 這些直方圖包含了以某個(gè)像素為中心的(2r+1)個(gè)像素信息(r為濾波窗口的半徑). 以每列直方圖第一個(gè)像素為中心的(2r+1)個(gè)像素信息為圖像初始狀態(tài), 更新處理時(shí), 這些直方圖隨滑塊移動(dòng)并完成更新.

2) 更新核直方圖. 隨著滑塊移動(dòng), 將每列直方圖減去最左側(cè)像素的直方圖信息, 如圖3(A)所示； 同時(shí)增加最右側(cè)像素的直方圖信息, 如圖3(B)所示, 每移動(dòng)一列都是O(1)操作.

圖3 更新核直方圖Fig.3 Update kernel histogram

3) 獲取最小值. 根據(jù)更新的核直方圖取最小值.

通過(guò)上述計(jì)算, 所有單像素操作都與窗口大小無(wú)關(guān), 即算法空間復(fù)雜度為O(1). 該算法偽代碼如下.

輸入： 圖像X, 大小為m×n, 濾波窗口半徑為r;

輸出： 與圖像X大小相同的圖像Y;

初始化核直方圖H及每列的直方圖h1,h2,…,hn

Fori=1 tomdo

Forj=1 tondo

RemoveXi-r-1,j+rfromhj+r

AddXi+r,j+rtohj+r

H←H+hj+r-hj+r-1

Yi,j← min(H)

End For

End For.

要將上述灰度圖像最小值濾波后擴(kuò)展適用于彩色圖像, 可分別對(duì)R,G,B 3個(gè)波段圖層最小值濾波, 然后比較3個(gè)波段圖層最小值濾波后的像素點(diǎn), 取3個(gè)像素值中的最小值作為暗通道值, 即實(shí)現(xiàn)了與濾波窗口無(wú)關(guān), 算法時(shí)間復(fù)雜度為O(1)的暗通道濾波.

使用該算法時(shí)間復(fù)雜度為O(1)的中值濾波器雖然對(duì)圖像邊緣保護(hù)效果并不突出, 但其是極簡(jiǎn)單且具有邊緣保護(hù)效果的濾波器. 實(shí)驗(yàn)表明, 在圖中云霧分布較均勻時(shí), 使用中值濾波器代替引導(dǎo)濾波效果較好; 當(dāng)云霧分布不均勻時(shí), 在暗通道先驗(yàn)去霧中使用中值濾波器的效果欠佳, 采用引導(dǎo)濾波. 基于暗通道先驗(yàn)的圖像去霧改進(jìn)算法流程如圖4所示.

圖4 暗通道先驗(yàn)去霧算法速度改進(jìn)流程Fig.4 Flow chart of speed improvement of dark channel prior fog removal algorithm

3 實(shí)驗(yàn)結(jié)果與分析

3.1 實(shí)驗(yàn)環(huán)境

當(dāng)圖像中云霧均勻分布時(shí), 用中值濾波代替引導(dǎo)濾波精細(xì)化大氣透射率圖層, 以提高精細(xì)化處理速度. 實(shí)驗(yàn)在湖北省恩施市航拍圖像中選取68幅云霧較均勻的圖像作為原始數(shù)據(jù)進(jìn)行驗(yàn)證. 圖片均為JPEG格式, 包含了R,G,B 3個(gè)顏色通道. 實(shí)驗(yàn)使用的計(jì)算機(jī)配置為： CoreTMi7-5500U @2.4 GHz處理器, 16 GB內(nèi)存, Windows10 64位操作系統(tǒng). 為使算法改進(jìn)前后的去霧效果偏差較小, 在不影響輸出圖像后期處理的基礎(chǔ)上盡量提高計(jì)算速度, 實(shí)驗(yàn)設(shè)定了統(tǒng)一參數(shù)批量處理全部原始有霧圖像.

3.2 改進(jìn)算法示例

實(shí)驗(yàn)使用無(wú)參考的客觀圖像評(píng)價(jià)方法, 將圖像信息熵、 邊緣強(qiáng)度、 方差作為評(píng)價(jià)指標(biāo). 以輸入如圖5所示的云霧均勻圖像為例, 文件大小7.3 MB, 使用中值濾波代替引導(dǎo)濾波在暗通道先驗(yàn)去霧算法中精細(xì)化處理大氣透射率圖層.

圖5 實(shí)驗(yàn)選取的有霧圖像Fig.5 Foggy image selected for experiment

圖6為分別使用中值濾波和引導(dǎo)濾波精細(xì)化的大氣透射率圖層, 其中光暈現(xiàn)象主要由計(jì)算暗通道時(shí)的塊狀處理結(jié)果導(dǎo)致. 由圖6可見(jiàn), 使用中值濾波精細(xì)化大氣透射率圖層效果不如引導(dǎo)濾波, 但中值濾波仍很好地模糊了因計(jì)算暗通道圖層而產(chǎn)生的塊狀處理結(jié)果, 并在一定程度保護(hù)了圖像邊緣信息.

圖6 使用中值濾波和引導(dǎo)濾波精細(xì)化后的透射率圖層Fig.6 Transmittance layer refined by median filtering and guided filtering

用圖6的兩幅大氣透射率圖層分別進(jìn)行暗通道先驗(yàn)去霧處理, 還原出的清晰圖像如圖7所示, 雖然使用中值濾波與引導(dǎo)濾波去霧后的兩幅清晰圖像色調(diào)存在一定差異, 但仍還原出了較清晰的圖像.

圖7 使用中值濾波和引導(dǎo)濾波去霧后的清晰圖像對(duì)比Fig.7 Comparison of clear images after defogging using median filtering and guided filtering

3.3 實(shí)驗(yàn)結(jié)果分析

實(shí)驗(yàn)選取68幅云霧較均勻的原始航拍有霧圖像(大小為5 180×5 400像素)進(jìn)行暗通道實(shí)驗(yàn), 其中部分原始有霧圖像如圖8所示. 將圖8原始圖像對(duì)應(yīng)的大氣透射率圖層使用中值濾波精細(xì)化, 得到如圖9所示優(yōu)化后的大氣透射率圖層.

圖8 原始航拍有霧圖像Fig.8 Original aerial foggy images

圖9 中值濾波處理后對(duì)應(yīng)的透射率圖層Fig.9 Corresponding transmittance layer after median filtering processing

用中值濾波將選取的68幅云霧均勻圖像經(jīng)過(guò)暗通道先驗(yàn)去霧后, 對(duì)實(shí)驗(yàn)輸出的全部清晰圖像進(jìn)行客觀圖像質(zhì)量評(píng)價(jià), 結(jié)果列于表1. 由表1可見(jiàn), 當(dāng)輸入圖像云霧均勻時(shí), 用中值濾波代替引導(dǎo)濾波對(duì)大氣透射率圖層精細(xì)化后的去霧效果有一定影響, 但仍能還原出清晰圖像. 當(dāng)去霧圖像應(yīng)用不受影響時(shí), 用中值濾波代替引導(dǎo)濾波精細(xì)化處理大氣透射率圖層顯著提高了暗通道先驗(yàn)去霧算法的運(yùn)算速度.

表1 精細(xì)化大氣透射率模塊改進(jìn)速度評(píng)價(jià)指標(biāo)Table 1 Evaluation indexes of refined atmospheric transmittance module for improving speed

綜上所述, 本文針對(duì)無(wú)人機(jī)航拍圖像中云霧情況多樣、 數(shù)據(jù)量大等特點(diǎn), 通過(guò)計(jì)算機(jī)編碼實(shí)驗(yàn), 提出了一種簡(jiǎn)單且具有邊緣保護(hù)效果、 與濾波窗口無(wú)關(guān)、 時(shí)間復(fù)雜度為O(1)的中值濾波器, 采用中值濾波代替引導(dǎo)濾波精細(xì)化處理大氣透射率圖層, 利用無(wú)參考的客觀圖像質(zhì)量評(píng)價(jià)方法量化評(píng)價(jià)算法改進(jìn)前后的去霧圖像效果. 結(jié)果表明, 當(dāng)輸入圖像中云霧較均勻時(shí), 中值濾波更易實(shí)現(xiàn), 且算法時(shí)間復(fù)雜度較低, 能節(jié)約去霧時(shí)間, 實(shí)現(xiàn)圖像平滑處理. 用中值濾波相比引導(dǎo)濾波去霧后的圖像效果有影響, 但仍可還原出清晰圖像, 在不影響應(yīng)用時(shí), 改進(jìn)算法的計(jì)算速度明顯提高. 由于中值濾波僅很好地模糊了計(jì)算暗通道圖層產(chǎn)生的塊狀處理結(jié)果, 并未真正加入圖像的梯度和紋理信息, 所以中值濾波對(duì)云霧分布不均的圖像處理效果不佳, 精細(xì)化處理仍需使用引導(dǎo)濾波.