2-Dimensional Perlin Noise Generator

I have recently finished updating my Perlin Noise Generator to work in two dimensions. This was actually a little bit more challenging that I expected, but here goes an explanation.

Perlin noise is an accumulation of many different wave functions, so the obvious thing to do is to separate each of these out and work on them one at a time. Each one needs a set of unique keypoints at regular intervals over the domain of the function. Because the noise generator is two dimensional, each keypoint needs to have an x coordinate, a y coordinate, and a magnitude. The amount of the keypoints should not be fixed, as the different functions that will need to be generated will have to have different frequencies and therefore different amounts of keypoints. Luckily, I wrote the one-dimensional keypoint generator with this in mind - now all I had to do to run this in two dimensions was to simply use tail recursion to call the one-dimensional keypoint generator many times (making sure to update the random generator each time, or all the random magnitudes will be the same) and append all the results to each other.

The next interesting part - figuring out which keypoints are the ones bounding a given (x, y) coordinate. The output of this function serves as input to an interpolation function, which uses the bounding box magnitudes to figure out what the magnitude of a given point would be. When I was thinking about how to do this, I was wondering what would happen if an (x, y) coordinate lied directly on a keypoint, or what if it doesn't lie in between a square of 4 but actually lies on a border of the square, in between exactly 2 keypoints. Based on this reasoning, I had this function return either 1, 2, or 4 keypoints; if the function returns 1 coordinate, the given point lies exactly on a keypoint, and the interpolation function can just use the magnitude of that point. If two points are returned, the function can use the linear interpolation that I already wrote in the last version of the noise generator. If four points are returned, then the full bicubic interpolation must be run.

That being said, my implementation of this function is not the best. In order to find if the given point lies on a keypoint, I filter the list of keypoints for one that lies on the given point. If this filtered list is empty, I proceed to filter the list of keypoints for all the ones that have only one coordinate equal to that coordinate in the given point (for the linear case). If this list isn't empty, I filter based on the greatest value that is less than the given point and the smalles value that is less than the given point, and return those. If none of those cases were true, I filter the list again twice to get all the points that are in the four quadrants around the given point. The top-left point is simply the last point in the top-left quadrant, and the bottom-right point is simply the first coordinate in the bottom-right quadrant (assuming the points come sorted, and filter doesn't rearrange the list). In order to get the top-right point and the bottom-left point, I had to do folds across the list to find one with minimum in one coordinate and maximum in the other. That makes for a ton of filters, each of which is O(n), and even 2 folds thrown into the mix. After I was done, running the program took a very long time.

Looking back, I feel that I should change the output of the function so that it always returns 4 points. If the given point lies on a point, I should return a specified square (say, to the bottom-right of the given point), and if a point lies on a line, I should return the square either below or to the right of the given point. This way, I can unify all the multiple cases into one case, limit the number of filters, and the interpolation function can deal with when a coordinate value = 0. Cubic and bicubic interpolation do this for us, so no additional code would have to be added.

The interpolation function itself is also a fairly interesting function. After reading the Wikipedia article on bicubic interpolation, it kindly provides a matrix that can be used. The use of this matrix is fairly simple - make a vector of all the values and partial derivatives as specified, multiply them by the given matrix, and use the output vector as coefficients into the simple double-sigma formula at the top of the page (Under "Bicubic spline interpolation"). Because Haskell is a very math-y language, implementing the double-sigma formula was extremely simple. After thinking for a second or two about how to do the matrix multiplication, I realized that Haskell's standard formulas were enough. Matrix multiplication is simply multiplying each entry i in row j with entry i in the vector. So, I am going to be doing some computation over ever row of the matrix - this looks like a map. The computation that I'm going to be doing is a dot product - which looks like a sum over a zip. All in all - the multiplication looks like this:

map (sum . (zipWith (*) vector)) matrix

Of course, this assumes that the matrix is in row-major form.

All in all, everything is mostly the same as the one-dimensional case, except with another dimension artificially added in. My code can be found here.

I did have a neat idea about making an image. All the two-dimensional noise images that I had seen have been grayscale, so that it looks kind of like fog. My idea was to simply run the algorithm 3 times, mapping the output into each of the red, green, and bue channels of the image, respectively. This would create a bright, more appealing image. Here is one such image, with 3 noise functions added together:

One dimensional Perlin Noise Generator

A few weeks ago, I became a little obsessed (oxymoron, I know) with Perlin Noise. It was originally developed by Ken Perlin as a way to create noise that has variation that varies. While that may sound like a tongue twister, it actually makes a lot of sense.

Simply taking a list of uniform random numbers (Note: this is *not* what Perlin noise is), and perhaps doing some Hermite Spline interpolation between each pair of adjacent points will make a landscape that is equally jaggy as it is flat. There will be no global variation in the random numbers - the "next" value in the list is entirely independent of all previous values.

Now, beautiful noise has the quality that the "next" value lies close to, but not exactly on the current value. This quality should be able to be observed at both a global and local environment. This means that, if you step back and squint, the noise should not be too jaggy or too flat. In addition, if you zoom in and look really close, you should still see small perturbations, and they should not be too grand but not too small that they don't exist either.

Perlin was able to achieve this using a pretty simple idea which appears to come from sound waves. Say you have a low (bass) tone and a high (treble) tone in a song. The low tone is going to have a low frequency but a large amplitude. The large magnitude is necessary because any wave with a low frequency is going to have small amount of energy carried with it, so the amplitude needs to be artificially increased to give the wave a decent amount of energy. The treble tone, however, has a high frequency and therefore does not need to have as high of an amplitude. Adding these two tones together (playing both tones at the same time) can sound like music, and will result in a wave that has small slight perturbations from a longer, more baseline (get the pun?) wave.

Now, simply adding two waves does not create the beautiful noise that we are in search of. We want to add together many, many waves, so we have slight perturbations on slight perturbations all the way down the line, and we get a nice, smooth, varying wave. The waves that we will add together, have a special property. Lets say we have one wave in our list of waves to add together. If the next wave we use has the property that it has twice the frequency and half the amplitude of the wave we already have, an interesting thing occurs. The rate of change of the two waves turns out to be exactly the same. Doubling the frequency halves the wavelength, and the amplitude is halved, so the entire wave is actually just being miniaturized (and different random values are being chosen, of course). So the resultant wave will then seem to vary smoothly, because the perturbations have the same average rate of change as all the other waves. In music, doubling the frequency is the idea of an octave, and playing multiple octaves together indeed sounds beautiful, the way we want our noise to sound.

I've been working on implementing a Perlin Noise generator in Haskell, and have come across a couple things to watch out for. The first is that, because Haskell is purely functional, you have to handle random numbers delicately. In an imperative language, every call to rand() has a side effect of changing the system's random seed. In Haskell, however, this side effect cannot occur. This means that the random function must take a random seed as an argument and return a random number and a split random seed. If you want to generate a new random number, you must then use the new random seed. Therefore, you have to be very careful about how your random seed evolves around the flow of your program. Any function that uses the random function must return a tuple with whatever it would return, and an updated random seed that the caller can continue with. It is quite a headache passing all these random seeds around, even though I can easily see why it is necessary.

I wanted to draw the output of the function with some image drawing library. I have used libGD many times in the past, in many languages (C, C++, PHP). After looking for a bit, I found that there is a Haskell module that provides bindings into libGD. I know that there is a Haskell package manager called Cabal, so I look for Cabal to install within Ubuntu's package manager. It's a no go; even a quick search on the Internet confirmed that Cabal is not currently in APT. Oh well, I'll just download the source myself and install it. This works like a charm on my Ubuntu box, and I look for documentation on how to use this module. There is basically none, and I spend at least a couple hours looking over the two or three useful pages on the Internet pertinent to the topic. Finally I scrap together a couple lines that does what I want it to do - simply draw some pixels on an image in a 'for' loop.

Then I try to run the code on my laptop. My laptop is running Snow Leopard, a 64-bit (mostly) OS. However, 64-bit support for GHC is not enabled on Snow Leopard, so I had to install the 32-bit version of GHC. However, now when I try to install the Haskell module, it must link to a 32-bit version of libGD. Okay, no problem, I go and run `port install libgd +universal` hoping that this builds the 32-bit version of libGD. It starts to, but then errors on some assembly file saying that some x86 assembly commands do not work on x86_64. Awesome. I can't download a binary from the libGD website, because I need the development libraries as that is what the Haskell module will hook into in the end. At this point I gave up; It already runs on my desktop.

But not for long! After a system update to Ubuntu, the darn thing prints out an error message about how libpthread.so has an invalid ELF header. At this point, I have had my fill trying to get this module to work, and I figure it would be much less of a headache just to write that specific part of the program in a language that is much widely accepted.

My idea is to have the Haskell script simply output the locations of all the pixels that should be turned on in the image to standard output. Then, my Perl script would gobble that up, draw the image with libGD, and output the image to standard output. If this works, I could draw the image by simply piping the two commands together. One Perl script later, this is working (And the Perl script is quite cool, may I say ... "while (<>)" is always a shocker when you see it =D Yay for default variables!)

While I am not quite done with this project (I want to implement 2-D and 3-D and 4-D versions of the same program), I have found a couple things that are wrong with Haskell. It's been around for at least 20 years, but the community just isn't there. There are a couple mailing lists that are quite large, and even some textbooks about it, but there was really no help I could find trying to use the little module that I was trying to use. Documentation is slim for these corner cases. Also, maybe I don't know the language well enough, but it seemed freaking impossible to do the simplest things with the image. Simply drawing pixels in a loop seemed to be much harder than it should be - I had to use a list comprehension to make a list of IO commands (one for drawing each pixel) and then make a "Do IO Commands" function that recursively runs through a list and does each command one-by-one. Dealing with IO commands, which is required for what I'm trying to do, just seems entirely awkward. Like I said earlier, though, maybe I just haven't learned how to use them elegantly yet.

Soon, I will post the same version of the program in many dimensions! look out for it!

My current code can be found here.