In my previous blog post, I talked about a personal project where I am generating heightmaps. I ended that post writing about how blurring has been important to me in generating the maps. Here is more information on the blurring.
First, I found out that rather then generating a lot of points, I could generate a lot less points and apply heavier blurring to get a satisfactory result much more rapidly.
For my blurring, I ended up using a Gaussian blur (sometimes called Gaussian filter). This type of filter is often used with height maps as it gives good natural looking results. Something very similar to erosion.
A Gaussian filter does this by giving more weight to the pixel being blurred rather than giving equal weight to the center pixel and each of it’s neighbors as is the case in the mean filter.
The problem I was facing was that my heightmaps needed to be blurred. Rather than create a blurring class that allows me to specify various degrees of blurring, work on images or even specify different types of blur filter (ie: gaussian, mean, etc.), I wrote the simplest thing that could solve my problem.
Applying YAGNI here has allowed me to complete this quickly and get on with the rest of the project. Years earlier I would have started a cool blurring library which I wouldn’t had the time to complete.
The Gaussian filter is a filter where the values of the kernel are calculated using the Gaussian function to produce values falling in a normal distribution. Like other filter (ie: the mean filter), the Gaussian filter works with a kernel which is a matrix.
At it’s simplest, a non-gaussian kernel could look something like this :
0.1, 0.1, 0.1,
0.1, 0.1, 0.1,
0.1, 0.1, 0.1,
We use the kernel to recalculate the value of each pixel of an image (or more simply of an array) by multiplying its actual value with the center position in the kernel and then multiplying the value above that pixel with the value above in the kernel and so on.
After we have we have multiplied all the nine values (for our 3 x 3 kernel), we add them up and store them in a single cell in a second array. We must do this for all points in our image or array.
When doing this, don’t forget to work with a secondary array so as to not modify the original. This is because you will be reading the original again when calculating the other pixels surrounding the first one you calculated.
For a Gaussian filter we would calculate the values of the kernel using a Gaussian distribution formula. I will not go into the details of the mathematical formula. You can find all the nitty gritty on the links at the beginning of the post or the references at the end.
Rather than create a 2d kernel, it is also possible to use a modified formula to create a 1d kernel. Running this 1d kernel twice, once horizontally and once vertically, produces a result equivalent to running the 2d kernel.
The difference is that running the 1d kernel twice is faster than running the 2d kernel once as it requires less operations. A neat trick.
The degree of blurring (how heavy to blur) is controlled by changing the number of standard deviations in our Gaussian distribution formula.
Rather than calculate a new kernel every time with differing values, another trick is to use the same kernel to blur our image many times in succession to obtain heavier blurring. If this is parallelized, this can also lead to speed optimizations.
For my implementation, I have used a precomputed 1d kernel I have found here and run it both horizontally and vertically.
Precomputed 1d kernel: [0.006, 0.061, 0.242, 0.383, 0.242, 0.061, 0.006].
If I need light blurring I will blur once and if I need heavier blurring (my most common scenario) I will blur twice in a row. Also I applied the blurring on an array of int and only created the image once the blurring was done which is much simpler than working with an image file.
Fast and simple to implement.
An issue I had to solve while implementing the blur was what to do when the filter/kernel is working with values on the edge of the array.
In these cases there are many possible solutions but I chose to use the value of the center pixel to fill in for the values of the non existing pixels.
The other choice would have been to wrap around the edge of the array, or to use a predefined value like 0 or 255. Using the center pixel gave pretty good results.
For the following code note that:
1- I am using a 1d array for my map to represent a 2d array which is unnecessary. If you prefer to use a 2d array, you should go with that.
2- The FilterOutOfBoundsSpecification which is in another file has been appended to the end of this code sample.
3- The code could be refactored as both compute_x_filtered_value and compute_y_filtered_value are almost identical and could be one function.
4- Even though the kernel is symmetrical, I store it likewise:
@filter_kernel = [[-3, 0.006], [-2, 0.061], [-1, 0.242], [0, 0.383], [1, 0.242], [2, 0.061], [3, 0.006]]
But I could save space by storing it this way:
@filter_kernel = [[0, 0.383], [1, 0.242], [2, 0.061], [3, 0.006]]
This would necessitate a small change to the code. I decided to keep the longer version of the kernel because I feel it’s more readable and easier to understand. Saving a few characters isn’t worth it if it complicates things.
5- For the sake of brevity I did not include the tests code.
Enough preamble here is the code:
require_relative './filter_out_of_bounds_specification' class GaussianFilter def initialize # the filter kernel value pairs. The first element in each pair represent # the distance from the central pixel and the second element the # weighted filter value @filter_kernel = [[-3, 0.006], [-2, 0.061], [-1, 0.242], [0, 0.383], [1, 0.242], [2, 0.061], [3, 0.006]] end # this function assumes the array represents a square 2d matrix def filter(array) @size = Math.sqrt(array.length) intermediate_filtered_array = Array.new(array.length, 0) filtered_array = Array.new(array.length, 0) # filter horizontally (0...@size).each do |y| # @size is used since we assume the array is square 2d matrix (0...@size).each do |x| intermediate_filtered_array[x + y * @size] = compute_x_filtered_value(array, x, y) end end # filter vertically (0...@size).each do |x| # @size is used since we assume the array is square 2d matrix (0...@size).each do |y| filtered_array[x + y * @size] = compute_y_filtered_value(intermediate_filtered_array, x, y) end end return filtered_array end private def compute_x_filtered_value(array, x, y) computed_value = 0.0 @filter_kernel.each do |filter_pair| if FilterOutOfBoundsSpecification.satisfied_by?(x, @size, filter_pair) offset = 0 else offset = filter_pair end computed_value += filter_pair * array[x + offset + y * @size] end return computed_value.round end def compute_y_filtered_value(array, x, y) computed_value = 0.0 @filter_kernel.each do |filter_pair| if FilterOutOfBoundsSpecification.satisfied_by?(y, @size, filter_pair) offset = 0 else offset = filter_pair end computed_value += filter_pair * array[x + (offset + y) * @size] end return computed_value.round end end class FilterOutOfBoundsSpecification def self.satisfied_by?(z, size_z, offset) return true if z + offset >= size_z return true if z + offset < 0 false end end
Here are the sources that I have used to write my Gaussian filter and this blog post, in order of relevance to my project and this article.
My post does not talk about the math and science side of the Gaussian distribution and about image processing as both subjects aren’t strengths of mine. So I encourage you to read the following links for more info on these subjects.