Bilateral Filtering¶

Important: Please read the installation page for details about how to install the toolboxes. $\newcommand{\dotp}[2]{\langle #1, #2 \rangle}$ $\newcommand{\enscond}[2]{\lbrace #1, #2 \rbrace}$ $\newcommand{\pd}[2]{ \frac{ \partial #1}{\partial #2} }$ $\newcommand{\umin}[1]{\underset{#1}{\min}\;}$ $\newcommand{\umax}[1]{\underset{#1}{\max}\;}$ $\newcommand{\umin}[1]{\underset{#1}{\min}\;}$ $\newcommand{\uargmin}[1]{\underset{#1}{argmin}\;}$ $\newcommand{\norm}[1]{\|#1\|}$ $\newcommand{\abs}[1]{\left|#1\right|}$ $\newcommand{\choice}[1]{ \left\{ \begin{array}{l} #1 \end{array} \right. }$ $\newcommand{\pa}[1]{\left(#1\right)}$ $\newcommand{\diag}[1]{{diag}\left( #1 \right)}$ $\newcommand{\qandq}{\quad\text{and}\quad}$ $\newcommand{\qwhereq}{\quad\text{where}\quad}$ $\newcommand{\qifq}{ \quad \text{if} \quad }$ $\newcommand{\qarrq}{ \quad \Longrightarrow \quad }$ $\newcommand{\ZZ}{\mathbb{Z}}$ $\newcommand{\CC}{\mathbb{C}}$ $\newcommand{\RR}{\mathbb{R}}$ $\newcommand{\EE}{\mathbb{E}}$ $\newcommand{\Zz}{\mathcal{Z}}$ $\newcommand{\Ww}{\mathcal{W}}$ $\newcommand{\Vv}{\mathcal{V}}$ $\newcommand{\Nn}{\mathcal{N}}$ $\newcommand{\NN}{\mathcal{N}}$ $\newcommand{\Hh}{\mathcal{H}}$ $\newcommand{\Bb}{\mathcal{B}}$ $\newcommand{\Ee}{\mathcal{E}}$ $\newcommand{\Cc}{\mathcal{C}}$ $\newcommand{\Gg}{\mathcal{G}}$ $\newcommand{\Ss}{\mathcal{S}}$ $\newcommand{\Pp}{\mathcal{P}}$ $\newcommand{\Ff}{\mathcal{F}}$ $\newcommand{\Xx}{\mathcal{X}}$ $\newcommand{\Mm}{\mathcal{M}}$ $\newcommand{\Ii}{\mathcal{I}}$ $\newcommand{\Dd}{\mathcal{D}}$ $\newcommand{\Ll}{\mathcal{L}}$ $\newcommand{\Tt}{\mathcal{T}}$ $\newcommand{\si}{\sigma}$ $\newcommand{\al}{\alpha}$ $\newcommand{\la}{\lambda}$ $\newcommand{\ga}{\gamma}$ $\newcommand{\Ga}{\Gamma}$ $\newcommand{\La}{\Lambda}$ $\newcommand{\si}{\sigma}$ $\newcommand{\Si}{\Sigma}$ $\newcommand{\be}{\beta}$ $\newcommand{\de}{\delta}$ $\newcommand{\De}{\Delta}$ $\newcommand{\phi}{\varphi}$ $\newcommand{\th}{\theta}$ $\newcommand{\om}{\omega}$ $\newcommand{\Om}{\Omega}$

This tour explores edge preserving filtering and in particular the bilateral filter, with applications to denoising and detail enhancement.

Gaussian Linear Filtering¶

The most basic filtering operation is the Gaussian filtering, that tends to blur edges.

Image size.

Load an image.

A Gaussian filter of variance $\si$ reads $$ \Gg_\si(f)(x) = \frac{1}{Z} \sum_y G_\si(x-y) f(y) \qwhereq Z = \sum_y G_\si(y), $$ and where the Gaussian kernel is defined as: $$ G_\si(x) = e^{-\frac{\norm{x}^2}{2\si^2}}. $$

A convolution can be computed either over the spacial domain in $O(N\si_x^2)$ operations or over the Fourier domain in $O(N \log(N))$ operations. Depending on the value of $ \si_x $, one should prefer either Fourier or spacial domain. For simplicity, we consider here the Fourier domain (and hence periodic boundary conditions).

Define a Gaussian function, centered at the top left corner (because it corresponds to the 0 point for the FFT).

Display a recentered example of Gaussian.

Define a shortcut to perform linear Gaussian filtering over the Fourier domain. This function is able to process in parallel a 3D block |F| by filtering each |F(:,:,i)|.

Example of filtering $ \Gg_\si(f_0) $.

Exercise 1

Display a filtering $\Gg_\si(f_0)$ with increasing size $\si$.

Bilateral Filter¶

The bilateral filter is a spacially varying filter that better preserves edges than the Gaussian filter.

It was first introduced in:

Carlo Tomasi and Roberto Manduchi, Bilateral Filtering for Gray and Color Images, Proceedings of the ICCV 1998

A very good set of ressources concerning the bilateral filter can be found on <http://people.csail.mit.edu/sparis/bf/ Sylvain Paris web page>. It includes a fast implementation, research and review papers.

Given a spacial width $\si_x$ and a value width $\si_v$, the filter opterates as: $$ \Bb_{\si_x,\si_v}(f)(x) = \frac{1}{Z_x} \sum_y G_{\si_x}( x-y ) G_{\si_v}(f(x)-f(y)) f(y) $$ where the normalizing constant is $$ Z_x = \sum_y G_{\si_x}( x-y ) G_{\si_v}(f(x)-f(y)).$$

At a given pixel $x$, it corresponds to an averaring with the data-dependant kernel $ G_{\si_x}( x-y ) G_{\si_v}(f(x)-f(y)) $.

Bilateral Filter by Stacking¶

Implementing the bilateral filter directly over the spacial domain requires $O( N \si_x^2 )$ operations where $N$ is the number of pixels.

A fast approximate implementation is proposed in:

Fast Bilateral Filtering for the Display of High-Dynamic-Range Images, Fredo Durand and Julie Dorsey, SIGGRAPH 2002.

It exploits the fact that for all pixels $x$ with the same value $v=f(x)$, the bilateral filter can be written as a ratio of convolution $$ \Bb_{\si_x,\si_v}(f)(x) = F_v(x) = \frac{ [G_{\si_x} \star ( f \cdot W_v )](x) }{ [G_{\si_x} \star W_v](x) } $$ where the weight map reads $$ W_v(x) = G_{\si_v}(v-f(x)). $$

Instead of computing all possible weight maps $W_v$ for all possible pixel values $v$, one considers a subset $\{v_i\}_{i=1}^p$ of $p$ values and computes the weights $ \{ W_{v_i} \}_i $.

Using these convolutions, one thus optains the maps $ \{ F_{v_i} \}_i $, that are combined to obtain an approximation of $ \Bb_{\si_x,\si_v}(f)$.

The computation time of the method is $ O(p N\log(N)) $ over the Fourier domain (the method implemented in this tour) and $ O(p N \si_x^2) $ over the spacial domain.

Value of $\sigma_x$:

Value of $\sigma_v$:

Number $p$ of stacks. The complexity of the algorithm is proportional to the number of stacks.

Function to build the weight stack $ \{ W_{v_i} \}_i $ for $v_i$ uniformly distributed in $[0,1]$.

Compute the weight stack map. Each |W(:,:,i)| is the map $W_{v_i}$ associated to the pixel value $ v_i $.

Exercise 2

Display several weights $ W_{v_i} $.

Shortcut to compute the bilateral stack $\{ F_{v_i} \}_i $.

Compute the bilateral stack $\{ F_{v_i} \}_i $.

Exercise 3

Display several stacks $F_{v_i}$.