鄒瑞波,廖海泳
(汕頭大學工學院,廣東 汕頭 515063)
高動態范圍成像(High Dynamic Range Imaging,簡稱HDRI或HDR),在計算機圖形學與電影攝影術中,是用來實現比普通數位圖像技術更大曝光動態范圍(即更大的明暗差別)的一組技術.高動態范圍成像的目的是要正確地記錄和重現真實場景中各區域的細節,特別是高亮區和暗區.隨著技術的發展,很多成像設備都可以達到很高的分辨率,例如手機拍照的分辨率已經可達2300萬像素,在這種情況下,高分辨率HDR圖像需要有一個有效和快速的算法.目前有各種各樣的HDR色調映射算法,然而很多算法不能同時兼顧保持細節信息和快速這兩個特點.例如,基于簡單變換函數的算法如文獻[1],雖然速度很快,但是由于壓縮亮度的同時也壓縮了細節的亮度,使得變換后邊緣細節信息不夠清晰.另一個方面,基于非線性變換的方法,如WLS算法[2]雖然能夠很好地重現細節信息,但是由于計算量極大,甚至不能實現,而不能應用在高分辨率HDR圖像.本文提出的色調映射算法首先利用移動最小二乘提取光照分量,也就是在移動窗口中利用最小二乘法提取光照分量.接著根據成像模型對光照分量進行壓縮,最后再恢復細節和色彩信息.實驗表明本文算法針對高分辨率HDR圖像能夠很好地保持和重現細節信息,而且速度比WLS算法快.
根據成像模型[3],圖像f(p)可以表示為:

其中p表示圖像像素的位置,i(p)為光照分量,r(p)稱為反射分量.即圖像在位置p的像素值由這兩者的乘積得到.光照分量i(p)的大小反應了在自然場景中照射到物體的光能的多少.根據我們的日常經驗可以知道,自然場景中的光照一般都是局部均勻,雖然有不同明暗的區別,但是從信號處理的角度看,所含的細節信息比較少.所以我們可以通過改變i(p)的值來壓縮圖像的亮度范圍,最終使得亮區和暗區的信息能夠同時較好地顯示出來.反射分量r(p)的幅值相比之下往往較小,但富含細節,能夠真實反映場景中的各種信息.因此這一分量的信息要盡量保存.
基于前面的理論分析,我們利用如下模型進行壓縮,重建LDR圖像:

其中f′(p)為壓縮重建后的LDR圖像.c為壓縮系數.根據冪函數的性質[4],我們可以知道當c>1時,將對暗區像素值進行壓縮,同時放大亮區像素值范圍;反過來,當c<1時,將對暗區像素值范圍進行放大,同時壓縮亮區像素值范圍,如圖1所示.通過大量的觀察我們發現絕大部分的HDR圖像像素值主要集中在灰度直方圖的左邊,而右邊的高像素值所占比例極少.圖2 Memorial的灰度直方圖就是一個經典的例子.綜上可知,為了達到增強暗區細節對比度的效果,c的設置一般要在小于1.

圖1 冪函數

圖2 Memorial灰度直方圖
我們先給出整個算法的主要步驟,再對關鍵步驟進行解釋.
算法步驟如圖3所示

圖3 算法框架
在第二部分中我們利用移動最小二乘方法提取光照分量.移動最小二乘方法在下一節中將進行詳細介紹.從圖2可以觀察到HDR圖像亮區的像素值比其他區大得多,而且亮區的像素在直方圖中所占比例非常小.根據這個特點,我們在第三部分中利用98%百分位數對高亮像素值進行截斷,這樣可以有效的縮短動態域范圍,進而達到壓縮的目的.
由于光照分量具有局部均勻的特點,如果采用傳統的平滑濾波器進行濾波會導致邊緣模糊,從而使得重建后的圖像中邊緣有光暈的問題.為了避免這個問題,我們可以采用各種邊緣保持濾波器進行提取,例如WLS濾波器[2],雙邊濾波器[5],TV濾波器[6]等.由于WLS濾波器具有很好的邊緣保持效果,下面主要介紹WLS濾波器的基本原理.

給定灰度圖像g(p),WLS濾波器通過求解來得到光照分量u(為了在下文避免符號容易混淆問題,這里使用u表示光照分量,而不是前面的i(p)).在式子(3)中第一項是數據保真項,用來保證u和g盡量接近,第二項是正則化項,用來控制所求u的平滑程度[7].是一個正則化參數, 越大則所求的光照分量u越平滑.而關于平滑權重的計算,我們參考了文獻[8]的方法,具體如公式4所示

其中I=log(g).ax,p(g),ay,p(g)分別表示在像素p位置求偏導數.α取值一般在[1.2,2]的范圍.ε是一個很小的值,例如,ε=10-4,甚至更小,主要用來避免零做除數.式子(3)是一個二次型,對其求導可得到一個線性方程組

其中矩陣

這里下標i,j表示圖像中的像素,N4(i)表示像素i的4-鄰域.u,g分別是光照分量和原圖像按列形成的列向量.系數矩陣A是一個大型的稀疏矩陣,其大小為N×N,對應圖像的像素個數為N.圖4是N=100時,系數矩陣A非零元素的分布模式.可以看到A是一個主對角帶狀的稀疏矩陣.值得一提的是,當圖像分辨率為1000×1000時,矩陣A的大小為106×106.而現在日常生活中圖像分辨率至少是這個級別,甚至更大,所以WLS濾波器求解一個大型稀疏線性方程組將會更加困難.
求解大型稀疏方程組有迭代法[9],直接法[10]兩大類.后面我們采用的是直接法求解.關于直接法,具體可參考方獻[11].

圖4 稀疏矩陣
基于WLS濾波器的色調映射方法具有很好的視覺效果,然而由于需要求解大型方程組,不適合高分辨率圖像.高分辨率圖像的像素個數一般可以達到107這個級別,如果直接對整個圖像進行求解線性方程組,計算量相當驚人,幾乎不可能實現.受平滑濾波器[12]和圖像分塊處理算法[13]啟發,對高分辨率圖像,我們采用平滑濾波器的框架,在移動窗口里面進行最小二乘濾波.具體如圖5所示,設置一個移動窗口,該窗口首先位于圖像的左上角,先在窗口里面利用式子(5)求出光照分量,再將窗口以步長h向右移動,接著求出新位置的光照分量,以此類推,一直到圖像的右邊邊界.然后窗口下移步長h,再從左到右移動濾波,重復之前步驟,直到窗口到達圖像的右下角.對于窗口移動過程中的重疊區域,我們采用加權平均的方式計算該重疊區域的像素值.

圖5 稀疏矩陣
為了驗證本文算法是否有效,我們將本文算的濾波效果與WLS濾波的效果進行比較,具體見圖6.其中圖6(a)測試圖片來自文獻[2].(b)(d)(f)是 WLS 濾波的效果,(c)(e)(g)是移動最小二乘濾波的效果,可以看出本文算法可以達到WLS濾波器[2]的效果,這說明了本文算法的有效性.


圖6 對比WLS和本文算法對于邊緣保持的效果
在移動最小二乘濾波器算法里面需要設置兩個參數,一個是移動窗口的大小,一個是移動窗口的步長(同時決定了相鄰窗口的重疊區域).下面討論這兩個參數的選擇對濾波效果和速度的影響.顯然,移動窗口的設置不能太小也不能太大,因為太小的話而導致重疊區域增多,從而影響整體的計算速度,相反,如果太大會增加在移動窗口中求解線性方程組的難度.移動窗口的步長會影響到相鄰窗口的重疊區域,例如,步長太大,則相鄰窗口的重疊區域太小,則可能會造成塊狀效應.反之,如果太小,則相鄰窗口的重疊區域較大,會增加一些不必要的計算量,從而影響計算速度.經過實驗我們發現,只要相鄰移動窗口有至少百分之十的重疊區域就可以達到WLS濾波器的效果.為了更好地測試窗口大小對計算速度的影響,我們采用不同大小的窗口對圖6(a)進行測試,在保證與WLS濾波器一樣效果和使用最長步長的情況下,表1列出了窗口大小和濾波計算時間的結果.可以看出當窗口大小為200×200時或以下時,計算速度較快.使用這樣大小窗口的另一個好處是:在這種情況下,求解線性方程組不需要占用大量的內存資源.

表1 不同窗口大小對計算時間的影響
在實驗中本文采用的硬件設備為:Core i5-7200U,@2.50Hz雙核處理器,Nvidia GeForce 920MX(2GB)顯卡,8GB內存.軟件環境為:Ubuntu16.04+MATLAB.由于篇幅所限,下面展示兩幅具有代表性的HDR圖像及其測試效果,其中圖6(a)Oxford Church HDR圖像的分辨率為1 013×1 200.圖7(a)Narrow Path HDR圖像的分辨率為8 192×4 096.圖7(b)到圖7(d)分別是本文算法,WLS 算法[2],雙邊濾波算法[5]處理后得到結果.可以看到,本文算法與WLS算法非常接近,而圖7(d)整體上細節沒有圖7(b)和圖7(c)清晰,特別是窗口玻璃上的細節,這是由于雙邊濾波器的邊緣保持效果沒有WLS濾波器好這一特點決定的,原因具體見文獻[2].同樣的,通過比較圖8(b)和圖8(c),我們可以觀察本文算法在細節方面明顯比雙邊濾波器清晰,特別是石頭后面遠處的山這些區域.需要說明的是,由于圖8(a)Narrow Path HDR圖像的分辨率過高,WLS濾波器需要求解一個超大規模的線性方程組而導致內存溢出,無法得到最終的結果.
最后,我們采用不同的HDR圖像進行測試,比較3個算法的處理時間,結果如表2所示.當圖像分辨率比較大時,本文算法的計算時間明顯比WLS算法少,而且能夠處理任意分辨率的HDR圖像,后者是WLS算法所做不到的.雖然雙邊濾波明顯比其他兩種算法快,但是從前面的實驗結果圖7,圖8可知基于雙邊濾波的算法整體上細節沒有比其他兩種算法清晰.

圖7 Oxford Church HDR圖像的算法效果比較

圖8 Narrow Path HDR圖像的算法效果比較
近年來,隨著成像設備分辨率的不斷提高,高分辨率HDR圖像受到越來越多人的青睞.WLS算法是一個很優秀的色調映射算法,然而對于高分辨率的HDR圖像,WLS算法由于需要求解超大規模的線性方程組而導致無法得到最終結果.本文提出的移動最小二乘濾波算法可以克服這一缺點,同時可以很好地保持圖像邊緣信息.實驗結果表明,本文中所提的色調映射算法能很好地保持圖像細節,能夠處理任意分辨率的HDR圖像.

表2 本文算法和WLS算法處理時間比較