Tone mapping: A quick survey. (Tone mapping 原理)

c1aa14506d4d0132e77d6ea4796472e0

假設已經有一個HDR影像,如何顯示在8-bit螢幕上哩? Tone mapping三巨頭,依照結果優排序! 三種都是local tonemapping,都是SIGGRAPH02耶! 那年真流行HDR。

  1. Gradient domain high dynamic range compression
  2. Fast Bilateral Filtering for the Display of High-Dynamic-Range Images
  3. Photographic Tone Reproduction for Digital Images

tm

2和3的方法皆屬於edge preserving,對比度較差;而方法1屬於gradient preserving,對比度較佳。當然這三種算法的結果也反映了演算法複雜程度...

 

 

在此之前假設已經有了HDR影像,參考合併多張LDR影像成HDR影像的方法

1. Gradient domain high dynamic range compression

1.1 Introduction

看到gradient domain,免不了要來個poisson reconstruction了(抱頭)。論文的想法很簡單,壓縮gradient的scale,並保持正負關係即可。那些巴拉巴拉human eye的假說,不管你信不信,反正我是信了!

注意這裡沒有從影像拆出luminance term(古早多半是分出低頻當luma),而是嘗試用gradient直接反推回影像,一氣呵成。所以我猜隔年(SIGGRAPH03)開始流行posisson image editing不是沒有原因的。

注意這邊只保護了gradient的正負關係,並沒有對影像的亮暗分布做保護,意思是調整後陰影區塊的平均亮度可能超過天空的亮度。但那又怎樣~~~拉拉拉?

1.2 Model

gradientcompress上圖順序是這樣看的

  1. 原訊號x
  2. 原訊號取log,H(x)=log(x)
  3. log後取gradient,H^{'}(x)
  4. H^{'}(x)壓縮一下,變成G(x)
    • 大gradient壓縮比較多,小gradient壓縮比較少
  5. G(x)重建新的log(x)
  6. exp(log(x))返回x

看起來很簡單,但困難點在步驟5─解PDE,白話點就是給你edge,求原影像。這可以解釋為一個non-blind deconvolution problem ( 已知point spread function是laplician ),論文用Full Multigrid其實是轉換成一個Ax=b的optimization problem。

  • Multigrid:要interpolation或restriction的流程已經寫死,如固定跑V cycle或W cycle。優點是不管有沒有convergence都會停止,適合硬體加速。
  • Full Multigrid:基本上就是依序走完所有scale的V cycle。但是有極大的redundance,改進方法如後述;要interpolation或restriction都取決於optimization的residual大小。運氣不好可能會無窮迴圈,需要小心設定判斷條件,但設定好可以保證是convergence的結果;不適合硬體加速。

1.3 Poisson equation derivative

Poisson equation描述了多維空間的gradient加總,同時考慮了所有維度的影響力

簡單認識一下Laplace operator,其中\nabla代表一次微分

\Delta f = \nabla^2 f = \nabla \cdot \nabla f

多維直角坐標空間\Delta f為各維度之加總

\Delta f = \sum_{i=1}^n \frac {\partial^2 f}{\partial x^2_i}

我們令影像 I=f ,在二維平面展開

\nabla^2{I} = \nabla_{x} \cdot \nabla_{x}{I}+\nabla_{y} \cdot \nabla_{y}{I}

我們在一次微分與二次微分中插入一個點對點的衰減函數\Phi,使得一次微分的scale變小,注意這邊衰減函數是可以共用的


\begin{align*}
\nabla^2{I}&=\nabla_{x} \cdot \Phi\cdot\nabla_{x}{I}+\nabla_{y} \cdot \Phi\cdot\nabla_{y}{I}\\
&=\operatorname{div}{G}
\end{align*}

到此就導出了著名的Poisson equation

1.4 Attenuation function

對gradient的衰減想法基本上是這樣的,對於很大的gradient給予較大的衰減,很小的gradient則給予較小的衰減,超小的gradient則給予些許強化,僅維持可以被目測細節的gradient即可。

衰減公式如下


\Phi\left(x,y\right)={\left(\frac{\left\|{\nabla}H\left(x,y\right)\right\|}{\alpha}\right)}^{\beta-1}

其中\left\|{\nabla}H\left(x,y\right)\right\|為gradient magnitude,\alpha可視為一個gradient threshold,比它大的則會讓除數大於1,反之則小於1,1-\beta次方得到一個如下圖的曲線進行gradient的調整。\alpha一般我們定義為影像平均gradient magnitude的0.1倍

Rendered by QuickLaTeX.com

1.5 Experiment

筆者試著用Wiener deconvolution解但效果太強烈,亮暗部嚴重反轉,猜想是因為該論文中的\operatorname{div}{G}是為Laplacian的近似,另外就是因為直接在一個scale下求解,而非pyramid導出,在frequency domain可能需要一些修正,還是 Ax=b pyramid optimization比較soft。可以參考這邊的code

左圖為Multigrid(optimize with conjugate gradient),右圖為Wiener deconvolution結果

deconvolution的方法

下面是conjugate gradient optimization的影像,可以見到由原本的edge慢慢填補回完整的影像。

2. Fast Bilateral Filtering for the Display of High-Dynamic-Range Images

這篇的改進很簡單,旨在利用bilaterial filter分離較佳的luminance與detail,藉此改善halo;方法則是跟以前幾篇大同小異。因為論文中要使用很大的kernel,其最大的貢獻則是利FFT減少了bilaterial filter的時間複雜度由O(n^2)降低為O(n\log{n})

由於bilaterial filter可以使得luminance在邊緣處較為sharp,因此halo現象可望減少。

NOTE: 2008年已有 Constant Time O(1) Bilateral Filtering

bilater

如上圖所示,將影像拆為large scale luminance和detail兩部分,利用gamma correction降低luminance scale,而後再合併detail即可。具體演算法如下:

  1. input intensity= 1/61*(R*20+G*40+B)
  2. r=R/(input intensity), g=G/input intensity, B=B/input intensity
  3. log(base)=Bilateral(log(input intensity))
  4. log(detail)=log(input intensity)-log(base)
  5. log (output intensity)=log(base)*compressionfactor+log(detail) - log_absolute_scale
  6. R output = r*10^(log(output intensity)), etc.

其中步驟1. 的luminance算法可以隨意改,步驟4. luminance的壓縮方法也可以再調整,gamma correction的對比度實在不好。

3. Photographic Tone Reproduction for Digital Images

3.1 Model

基本上這個方法只能壓縮亮部,對暗部的對比提升沒有太大的幫助,但優點就是該方法沒有太複雜的數學,適合應用。如下式

\begin{equation}\label{eq:3.1}
L_{d}=\frac{L}{1+L}
\end{equation}

其中L_{d}為輸出的intensity,值域為[0, 1]L為輸入的HDR pixel。當L很大時,分母的1幾乎可以忽略,L_{d}輸出的值將會很接近1。反之當L很小時,L_{d}則接近直接輸出。如下圖

Rendered by QuickLaTeX.com

最後再乘上255,如此一來就可以保證輸出的intensity一定是在顯示器範圍。

3.2 Modified model for local tonemapping

\eqref{eq:3.1}是一個global tonemapping,分母的L相當於一個控制scale的常數,為了加強亮部細節,論文中將分母的L以區域亮度V取代。

\begin{equation}\label{eq:3.2}
L_{d}=\frac{L}{1+V}
\end{equation}

3.3 local luminance V

尋找區域亮度V以強化亮部壓縮的細節就是該論文重點。論文中以DoG (difference of Gaussian)實作。VL相差較大時則會保留細節。

這裡用到Blob detection的技巧,利用DoG或LoG (Laplician of Gaussian)來偵測一個合適的feature scale (在SIFT中最廣為人知),在這裡是尋找一個以L為中心的最大區塊使得內部的亮度最接近,並且該區塊沒有顯著的gradient變化。

我不喜歡論文前後不一致的的符號,所以稍作修改。

定義N個不同Gaussian kernel G_{i}模糊過的的影像B_{i},Gaussian的大小由\sigma_{i}控制。N在論文中為8,\sigma_{i}倍率在論文中為1.6倍,kernel大小由1至43 pixels,模糊程度由小到大排序。

B_{i}=L{\otimes}G_{i},{\quad}i=1{\dots}N

然後點對點相減,計算N-1個normalized DoG影像,其中c為跟\sigma_{j}有關的常數,c基本上是一個\sigma_{j}^2的倒數,用來防止B_{j}(x,y)接近 0 的時候會出問題

D_{j}(x,y)=\left|\frac{B_{j}(x,y)-B_{j+1}(x,y)}{c+B_{j}(x,y)}\right|,{\quad}j=1{\dots}N-1

對每個pixel (x,y)找一個最大的image index k(x,y)滿足下式

\max_{k(x,y){\in}j}{\left\{{k(x,y)}|{D_{k(x,y)}(x,y)<\varepsilon}\right\}}

最後得到區域亮度的影像V,其中每個pixel的值V(x,y)代表以該位置為中心,範圍最大的區域亮度,這個範圍中沒有明顯的gradient改變,感覺有點類似MSER

V(x,y)=B_{k(x,y)}(x,y)

3.4 Paper revision

這篇論文雖然方法簡單,但是要正確理解式子卻很難;照著式子寫程式一定錯的,因為該論文中有諸多刊誤和undocumented parameters (例如\alpha_{i}),要心領神會作者的表達需要一些想像力,或是作者給各位的彩蛋...

  1. 在Initial luminance mapping的式(1),\frac{1}{N}應該放在括號裡面,才符合文字敘述的log-average

    \bar{L}_{w}=exp\left(\frac{1}{N}\sum_{x,y}{log(\delta+L_{w}(x,y))}\right)

4 Replies to “Tone mapping: A quick survey. (Tone mapping 原理)”

  1. 其最大的貢獻則是利FFT減少了bilaterial filter的時間複雜度由O(n)降低為O(nlogn)。 => 其最大的貢獻則是利FFT減少了bilaterial filter的時間複雜度由O(n^2)降低為O(nlogn)。

Leave a Reply

Your email address will not be published. Required fields are marked *