Voronoi Stippling

I'm trying to write my own software to convert normal raster images to a corresponding stippled version.  If you're unfamiliar with the term, stippling is a style of artwork that that uses are large number of small dots to represent a scene.  Areas with a high density of dots appear darker, while areas without dots are white.  Normally I like to have reasonably good understanding of a subject before I do a blog post about it, but in this case I'll give a quick outline of what I'm trying to achieve and then fill in the details in a later post.

So, how do you get an image like this,

Stippled Image

From something like this.

Original Image

To start with I'm trying to replicate the process in a paper by Adrian Secord.

Weighted Voronoi Stippling
Adrian Secord
2nd International Symposium on Non-Photorealistic Animation and Rendering (NPAR 2002)

The first step is to generate a number of points equal to the desired number of stippling points in the final image.  They can be randomly distributed over the image, or ideally via a rough approximation, close to the locations of the stippling dots in the output image.  As the process I'm about to describe is iterative, the closer the initial approximation is the final output, the faster the algorithm will converge.

After the points have been created, their Voronoi diagram is generated.  A Voronoi diagram produces a region for each input point that identifies all the locations that are closest to that dot than any other.  There are some tricks to making sure that the regions around the boundary are properly created but in general it's a simple call to scipy.spatial.Voronoi

That's great, but the process so far hasn't taken the input image into account(unless you use it to generate the initial point locations).  Now that a number of points and their corresponding regions are created, the weighted centroid of the area is calculated.  The centroid then becomes the location of the new point.  This process of creating Voronoi diagrams and calculating centroids is repeated until some convergence criteria is met and is better known as Lloyd's Algorithm.

In the image below, the black dots are the initial points and the white dots are the centroids of the Voronoi regions after an iteration.

Voronoi Diagram with centroids

Overlaying this diagram on the original image helps to illustrate the process a little better.  Obviously the image below is sparsely stippled and the object isn't recognisable with only 64 points.

Voronoi Iteration

I've generated the next diagram with I think 256 points to explain a subtle point.  The Voronoi region for each point is a convex polygon.  Therefore its centroid will always be located within itself.  As the number of points increases the size of the regions decreases.  This means that the process will take a long time to converge as the points can only move a small amount during each iteration.   That's why you need to place the initial points as close to their final location as possible.  Or do you?

I only discovered this problem after programming my solution and needed a quick way to overcome it.  If a large number of points take a long time to converge then a small number should converge quickly.  So I started with only two points.  After a few iterations, each point splits, like a cell undergoing Mitosis.  Each point splits into two new points diametrically opposed on a unit circle centred on the old point.  The direction of the split is randomly chosen.  (I have a hunch that continually splitting in the same direction will cause artifacts that will delay convergence)  These new points are then iterated on for a while and then split again.  The key is to get the points close to where they need to be early and then converge and add the detail later.

Dense Voronoi diagram

I have some ideas I want to try to make the image better.  There are obviously some linearity problems with the output.  White parts aren't as light as they are should be.  My code is a mess too.  It needs a clean up.  I'm having fun though.