Laplacian of gaussian filter use

As you've probably figured out by now from the other answers and links, LoG filter detects edges and lines in the image. What is still missing is an explanation of what σ is.

σ is the scale of the filter. Is a one-pixel-wide line a line or noise? Is a line 6 pixels wide a line or an object with two distinct parallel edges? Is a gradient that changes from black to white across 6 or 8 pixels an edge or just a gradient? It's something you have to decide, and the value of σ reflects your decision — the larger σ is the wider are the lines, the smoother the edges, and more noise is ignored.

Do not get confused between the scale of the filter (σ) and the size of the discrete approximation (usually called stencil). In Paul's link σ=1.4 and the stencil size is 9. While it is usually reasonable to use stencil size of 4σ to 6σ, these two quantities are quite independent. A larger stencil provides better approximation of the filter, but in most cases you don't need a very good approximation.

Just for visualization purposes, here is a simple Matlab 3D colored plot of the Laplacian of Gaussian (Mexican Hat) wavelet. You can change the sigma(σ) parameter and see its effect on the shape of the graph:

sigmaSq = 0.5 % Square of σ parameter
[x y] = meshgrid(linspace(-3,3), linspace(-3,3)); 
z = (-1/(pi*(sigmaSq^2))) .* (1-((x.^2+y.^2)/(2*sigmaSq))) .*exp(-(x.^2+y.^2)/(2*sigmaSq)); 
surf(x,y,z)

You could also compare the effects of the sigma parameter on the Mexican Hat doing the following:

t = -5:0.01:5;
sigma = 0.5;

mexhat05 = exp(-t.*t/(2*sigma*sigma)) * 2 .*(t.*t/(sigma*sigma) - 1) / (pi^(1/4)*sqrt(3*sigma));

sigma = 1;
mexhat1 = exp(-t.*t/(2*sigma*sigma)) * 2 .*(t.*t/(sigma*sigma) - 1) / (pi^(1/4)*sqrt(3*sigma));

sigma = 2;
mexhat2 = exp(-t.*t/(2*sigma*sigma)) * 2 .*(t.*t/(sigma*sigma) - 1) / (pi^(1/4)*sqrt(3*sigma));

plot(t, mexhat05, 'r', ...
     t, mexhat1, 'b', ...
     t, mexhat2, 'g');

Or simply use the Wavelet toolbox provided by Matlab as follows:

lb = -5; ub = 5; n = 1000; 
[psi,x] = mexihat(lb,ub,n); 
plot(x,psi), title('Mexican hat wavelet')

I found this useful when implementing this for edge detection in computer vision. Although not the exact answer, hope this helps.

This was something that confused me too, and it wasn't until I had to do the same as you for a uni project that I understood what you were supposed to do with the formula!

You can use this formula to generate a discrete LoG filter. If you write a bit of code to implement that formula, you can then to generate a filter for use in image convolution. To generate, say a 5x5 template, simply call the code with x and y ranging from -2 to +2.

This will generate the values to use in a LoG template. If you graph the values this produces you should see the "mexican hat" shape typical of this filter, like so:

_{(source: ed.ac.uk)}

You can fine tune the template by changing how wide it is (the size) and the sigma value (how broad the peak is). The wider and broader the template the less affected by noise the result will be because it will operate over a wider area.

Once you have the filter, you can apply it to the image by convolving the template with the image. If you've not done this before, check out these few tutorials. java applet tutorials more mathsy.

Essentially, at each pixel location, you "place" your convolution template, centred at that pixel. You then multiply the surrounding pixel values by the corresponding "pixel" in the template and add up the result. This is then the new pixel value at that location (typically you also have to normalise (scale) the output to bring it back into the correct value range).

The code below gives a rough idea of how you might implement this. Please forgive any mistakes / typos etc. as it hasn't been tested.

I hope this helps.

private float LoG(float x, float y, float sigma)
{
    // implement formula here
    return (1 / (Math.PI * sigma*sigma*sigma*sigma)) * //etc etc - also, can't remember the code for "to the power of" off hand
}

private void GenerateTemplate(int templateSize, float sigma)
{
    // Make sure it's an odd number for convenience
    if(templateSize % 2 == 1)
    {
        // Create the data array
        float[][] template = new float[templateSize][templatesize];

        // Work out the "min and max" values. Log is centered around 0, 0
        // so, for a size 5 template (say) we want to get the values from
        // -2 to +2, ie: -2, -1, 0, +1, +2 and feed those into the formula.
        int min = Math.Ceil(-templateSize / 2) - 1;
        int max = Math.Floor(templateSize / 2) + 1;

        // We also need a count to index into the data array...
        int xCount = 0;
        int yCount = 0;

        for(int x = min; x <= max; ++x)
        {
            for(int y = min; y <= max; ++y)
            {
                // Get the LoG value for this (x,y) pair
                template[xCount][yCount] = LoG(x, y, sigma);
                ++yCount;
            }
            ++xCount;
        }
    }
}

It appears to be a continuous circular filter whose radius is sqrt(2) * sigma. If you want to implement this for image processing you'll need to approximate it.

There's an example for sigma = 1.4 here: http://homepages.inf.ed.ac.uk/rbf/HIPR2/log.htm