Heat equation $\frac{\partial u}{\partial t}=\Delta u$ is extremely useful for smoothing images,
where $u(x, 0)=f(x)$ i.e. $u$ is initialized to the given original image $f$, the image to be smoothed, we also impose Nuemann boundary conditions. Smoothing an image removes noise in image $f$. The problem with heat equation is though that it can not distinguish between noise and useful features like edges. Plus you never know when to stop! It just keeps on smoothing.
The boundary condition that we could use is Neumann boundary condition. i.e. the normal derivatives at the boundary is zero. The way to maintain this is to reflect the boundary, run the heat equation in the interior of the image, then shrink it back to the original size and repeat.
The Laplacian operator $\Delta u=u_{xx}+u_{yy}$ is easy to discretize:
$\Delta u \equiv Lu=(u_{i, j-1}+u_{i-1, j}+u_{i, j+1}+u_{i+1, j}-4u_{i, j})/h^2.$
It's a five point operator,
$$
Lu=\frac{1}{h^2}\left[\begin{array}{ccc}0 & 1 & 0 \\1 & -4\,\,\,\,\,\, & 1 \\0 & 1 & 0\end{array}\right]
$$
Before we perform the heat equation on an image $f$, we should see how the image really looks like.
>> f=imread('eyes.jpg');
>> figure; imshow(f);
This is what we get...
Let us try to see view the image as a graph of a function. This image being an RGB image we have to change it to grayscale before viewing it.
>> figure; mesh(double(rgb2gray(f)));
This command gives the following output.
In the above image, the magnitude of the grayscale is given by with blue indicating small magnitude and red color indicating large magnitide. This is a very "spiky" image. If we run heat equation on this image we should get a smoother image.
Here is the code I wrote for the heat equation. By now I assume that you have been convinced that the for loops are bad, so one can write the code heat_2fast without for loops.
>> u=heat_2fast(f, 1, 0.25, 1); figure; imshow(u/255);
Let us try to see view the image as a graph of a function:
Simply put heat equation diffuses the function f. As a matter of fact it diffuses the function f iso-tropically. i.e. equally in all directions. Let us see what happens if we run the heat equation for t =1000. This is what I got as a result:
This really flattened things out, huh? If we see this as a graph of a function we get the following.
This is what happens if you run the heat equation for t=1000. If you keep it running for more time eventually things become really flat like Kansas.
where $u(x, 0)=f(x)$ i.e. $u$ is initialized to the given original image $f$, the image to be smoothed, we also impose Nuemann boundary conditions. Smoothing an image removes noise in image $f$. The problem with heat equation is though that it can not distinguish between noise and useful features like edges. Plus you never know when to stop! It just keeps on smoothing.
The boundary condition that we could use is Neumann boundary condition. i.e. the normal derivatives at the boundary is zero. The way to maintain this is to reflect the boundary, run the heat equation in the interior of the image, then shrink it back to the original size and repeat.
The Laplacian operator $\Delta u=u_{xx}+u_{yy}$ is easy to discretize:
$\Delta u \equiv Lu=(u_{i, j-1}+u_{i-1, j}+u_{i, j+1}+u_{i+1, j}-4u_{i, j})/h^2.$
It's a five point operator,
$$
Lu=\frac{1}{h^2}\left[\begin{array}{ccc}0 & 1 & 0 \\1 & -4\,\,\,\,\,\, & 1 \\0 & 1 & 0\end{array}\right]
$$
Before we perform the heat equation on an image $f$, we should see how the image really looks like.
>> f=imread('eyes.jpg');
>> figure; imshow(f);
This is what we get...
 |
| Image f |
Let us try to see view the image as a graph of a function. This image being an RGB image we have to change it to grayscale before viewing it.
>> figure; mesh(double(rgb2gray(f)));
This command gives the following output.
 |
| Image f viewed as a graph of a function |
In the above image, the magnitude of the grayscale is given by with blue indicating small magnitude and red color indicating large magnitide. This is a very "spiky" image. If we run heat equation on this image we should get a smoother image.
Here is the code I wrote for the heat equation. By now I assume that you have been convinced that the for loops are bad, so one can write the code heat_2fast without for loops.
If you write the following in the command window:
>> u=heat_2fast(f, 1, 0.25, 1); figure; imshow(u/255);
 |
| Output of the heat equation at t=1. |
Let us try to see view the image as a graph of a function:
>> v=rgb2gray(u/255);
>> mesh(v);
 |
| Result of heat equation, it makes the image smooth ! |
The reason why heat equation is called a diffusion equation
Simply put heat equation diffuses the function f. As a matter of fact it diffuses the function f iso-tropically. i.e. equally in all directions. Let us see what happens if we run the heat equation for t =1000. This is what I got as a result:
 |
| The result of heat equation after a long time |
This really flattened things out, huh? If we see this as a graph of a function we get the following.
 |
| Heat equation FTW! |
This is what happens if you run the heat equation for t=1000. If you keep it running for more time eventually things become really flat like Kansas.
