Arkanis Development

A few days ago I needed a C function to blur an image. I had an old implementation of a Gaussian blur filter lying around so I reused that. It worked rather well and it's just one C function so I decided to clean it up and package it into a single header file library: iir_gauss_blur.h. This is what it looks like:

int width = 0, height = 0, components = 1;
uint8_t* image = stbi_load("foo.png", &width, &height, &components, 0);

float sigma = 10;
iir_gauss_blur(width, height, components, image, sigma);

Here stb_image.h is used to load the image before it is blurred. iir_gauss_blur.h compiles as C and C++ so you can just plug it in and use it (tested it on GCC only though).

The function is an implementation of the paper "Recursive implementation of the Gaussian filter" by Ian T. Young and Lucas J. van Vliet. It has nothing to do with recursive function calls, instead it's a special way to construct a filter. Other (convolution based) gauss filters apply a kernel for each pixel and the kernel grows as the blur strength (sigma) gets larger. Meaning their performance degrades the more blurry you want your image to be.

Instead the algorithm in the paper sweeps across the image four times: From left to right, right to left, top to bottom and bottom to top. This is enough to propagate the blur across the entire image, no matter how large or small the blur strength (sigma) is. Thanks to that the performance always stays the same. You can have some ridiculously large sigmas without any performance hit.

The meaning of life… eh, sigma (σ)

Sigma (σ) is a number (at least 0.5) that defines the blur strength:

It seems the usual approach is to start out with some sigma value and then adjust it until you get the blurriness you want.

There's also the concept of a "blur radius" but I have no idea what exactly it's supposed to mean. I suspect it's the kernel size for convolution based algorithms. There seem to be several rules of thumb out there, e.g. radius = 2 * sigma. So if you want to have a blur radius of 10px you can use sigma = (1.0 / 2.0) * radius to get the sigma for it (5.0). I also found other interesting equations here and here (again for kernel sizes I guess).

All this got me thinking about what a useful "blur radius" would actually mean for me. Mentally I think that the filter blurs each pixel and I want to configure the radius by which a pixel is blurred. In my current project I have two difference scenarios where that's important for me:

• How far a pixel is blurred visually. That's the distance from the pixel where it can still leave a visible impact (cause a visible change). For "visible change" I used the somewhat arbitrary limit of a change of 16 in the brightness of a pixel (0 to 255). I used that since that's the smallest change from black I could make out on my display.
• How far we have to be away from a pixel until it can't change the color values at all. That's the distance at which it's impact isn't strong enough to change an 8 bit pixel value (0 to 255) by even 1. It's useful for me since that's the distance I need when creating padding around an image so I can be sure that the border of the blurred image stays black.

I tried to solve the equation of the normal distribution to calculate the proper sigma for a given radius and threshold (e.g. 16/256 or 1/256 as above). But after a bit of tinkering (and a lot of math confusion) I gave up and went the empirical way: I made a test image and blurred it with different sigmas. Then used GIMP to measured how far the white seeped into the black region (again for the thresholds of 16 and 1). This is what I got for the threshold of 16:

I did the measurements one image (and sigma) at a time. Combining many into the image above makes the lower gray values a bit more visible than they're otherwise. Anyway, I did the same for the threshold of 1 and plugged it into Wolfram Alpha to approximate the values with a linear equation. Again, after some fiddling around this is what I got:

sigma = (1.0 / 1.42) * radius16
sigma = (1.0 / 3.66) * radius1

It seems to work more or less well for a sigma up to 100. At least that's good enough for my use cases.

But given how crude these approximations are take them with a grain of salt. Or even better a whole tee spoon of salt. They work for me but might not work for anything you have in mind. That's why I decided to leave these calculations out of the iir_gauss_blur() function.

Just in case they work for you and you also want to use them via the command line I created a corresponding CLI tool.

By the way: The algorithm looks like it could be implemented quite nicely on a GPU with OpenCL or compute shaders. I assume for small sigmas a convolution based GLSL implementation is probably faster. It can use the bilinear interpolation of the hardware and is probably more local and texture cache friendly. But for larger sigmas the "recursive" algorithm will outperform it, even if it needs 4 passes over the data. Also its performance doesn't degrade with the blur strength so that might be the killer argument in some cases. Anyway, I'm very happy with the (single-threadded) CPU performance so I won't do that any time soon.

I hope this stuff is useful to someone or was at least interesting. :)