Linear Diffusion

Page last edited 4,156 days ago
From Agung Alfiansyah
Jump to navigation Jump to search

Linear Diffusion Flows

This tours studies linear diffusion PDEs, a.k.a. the heat equation. A good reference for diffusion flows in image processing is [Weickert].

Contents

Installing toolboxes and setting up the path. Heat Diffusion Explicit Solution using Convolution Bibliography Installing toolboxes and setting up the path.

You need to download the following files: signal toolbox and general toolbox.

You need to unzip these toolboxes in your working directory, so that you have toolbox_signal and toolbox_general in your directory.

For Scilab user: you must replace the Matlab comment '%' by its Scilab counterpart '//'.

Recommandation: You should create a text file named for instance numericaltour.sce (in Scilab) or numericaltour.m (in Matlab) to write all the Scilab/Matlab command you want to execute. Then, simply run exec('numericaltour.sce'); (in Scilab) or numericaltour; (in Matlab) to run the commands.

Execute this line only if you are using Matlab.

getd = @(p)path(p,path); % scilab users must *not* execute this Then you can add the toolboxes to the path.

getd('toolbox_signal/'); getd('toolbox_general/'); Heat Diffusion

The heat equation reads for a function and where (the solution at initial time ) is given.

The Laplacian operator reads Δ The flow is discretized in space by considering a discrete image of pixels.

n = 256; Load an image , that will be used to initialize the flow at time .

name = 'hibiscus'; f0 = load_image(name,n); f0 = rescale( sum(f0,3) ); Display it.

clf; imageplot(f0);

The flow is discretized in time using an explicit time-stepping Δ We use finite difference Laplacian Δ where we assume periodic boundary conditions, and where is the spacial step size.

h = 1/n; delta = @(f)1/h^2 * div(grad(f)); The step size should satisfy for the discretized flow to be stable.

The discrete solution converges to the continuous solution at time if both and under the condition .

Select a small enough step size.

tau = .5 * h^2/4; Final time.

T = 1e-3; Number of iterations.

niter = ceil(T/tau); Initialize the diffusion at time .

f = f0; One step of discrete diffusion.

f = f + tau * delta(f); Exercice 1: (the solution is exo1.m) Compute the solution to the heat equation.

exo1;

Explicit Solution using Convolution

The solution to the heat equation can be computed using a convolution where denotes the convolution of continuous functions and is a Gaussian kernel of width

One can thus approximate the solution using a discrete convolution. Convolutions can be computed in operations using the FFT, since

cconv = @(f,h)real(ifft2(fft2(f).*fft2(h))); Define a discrete Gaussian blurring kernel of width .

t = [0:n/2 -n/2+1:-1]; [X2,X1] = meshgrid(t,t); normalize = @(h)h/sum(h(:)); h = @(t)normalize( exp( -(X1.^2+X2.^2)/(4*t) ) ); Define blurring operator.

heat = @(f, t)cconv(f,h(t)); Example of blurring.

clf; imageplot(heat(f0,2));

Exercice 2: (the solution is exo2.m) Display the heat convolution for increasing values of .

exo2;

Bibliography


[Weickert98] Joachim Weickert, Anisotropic Diffusion in Image Processing, ECMI Series, Teubner-Verlag, Stuttgart, Germany, 1998.

Copyright (c) 2010 Gabriel Peyre