Islands Are Just Mountains Up To Their Necks in Ocean

I've posted previously about creating islands.  In keeping with the mountain theme on Reddit, I'm going to talk about how to create archipelagos out of mountains.

Some islands -- like the Hawaiian islands -- are formed by volcanoes rising from isolated spots on the sea floor.  Other islands -- like the Caribbean islands -- are formed by the same sort of plate tectonics that create mountain ranges.  It's just that the action happens where the ocean is deep enough to cover up most of the mountains.  Only the tops of the "mountain range" rise out of the sea, forming a chain of islands along the edge of the tectonic plate.  This suggests creating an island chain by using the same process used to create mountain ranges but putting them in an ocean area and sinking them so that only the tops are visible.

While this would work to create free-standing ocean islands like the Caribbean islands, my maps don't have a lot of open ocean, so I'm going to use this technique to create island chains that extend off the end of a peninsula, like the Florida Keys:

Or the South Shetland Islands off the coast of Antarctica:

This actually happens "naturally" on some of the maps that Dragons Abound creates, where a mountain range runs into the sea.  But because of the way Dragons Abound creates ocean areas, it doesn't happen often, and I'd like a way to purposely create archipelagos.

The first problem is to create the proper shape of islands -- a kind of wedge-shaped area of islands stretching out from a peninsula.  So if I have a bit of map that has a little peninsula like this:

I want to add some islands in roughly this area:

If the islands start dense and wide near the shore and get narrower and lower as they stretch out to see, it will hopefully appear as if the same ridge that created the peninsula continues out to sea, sinking lower and becoming islands.

As it happens, the same tent mask that I used to create mountain ranges will serve for this purpose as well. 

The end of the tent mask creates exactly the wedge shape that I need.  I just need to set the length of the tent to zero and just use one end.  To test this out, I'll try filling the mask with a solid chunk of land. (I'm placing this off the edge of a simple circular land shape to get a notion of how it would sit in relation to a peninsula.)

That looks pretty good.  Now I need to generate mountains in the mask rather than just flat land.

Here I've used noise parameters that will give me a lot of small pointy mountains -- three octaves of ridged multi-fractal noise.  I chose this because I want to get a lot of small islands, but you can play around with different noise formulae to find a setting that you like.

Now I need to sink the mountains into the ocean so that only the tops are showing.  With a little playing around with heights, I got this:

Which already looks pretty good.  I control how far to sink the mountains by setting a percentage of area within the mask that should be islands, and then I lower all the mountains until only that percentage is above water.  (I'm showing more islands in these images than I think looks best on a map because I want the shape of the chain to be evident.  In practice, I find about 25% land looks good.)

For a further refinement, you can perturb the shape of the mask with noise to get something that is less regular.  Here is the basic mask:
And here are some perturbed versions: 

I'm using fairly minor perturbations here.  It's important that the island chain point away from peninsula to be convincing, so I don't want to lose that basic arrow shape.  Here is an example of an island chain created with a perturbed mask:

Obviously it would be easy to use other mask shapes here, such as a long narrow arc to create something like the Caribbean islands.

Now that I can create an island chain, the next step is to identify candidate peninsulas on the map and figure out how to properly align the mask with the peninsula.  Sadly, the research field of peninsula identification seems to be moribund, so I was forced to invent my own algorithm.

A peninsula is a point on the coast line where the land bulges out into the sea.  To detect these spots, I slide 2 points (Beginning and End) along the coast line, keeping them a fixed distance along the coast apart.  At each step, I check whether the midpoint of the line between Beginning and End is over land.  For example:

Here the green point is on the midline between Beginning and End and is over land.  So I know that the coast between Beginning and End forms a peninsula.  If there is water, I have a bay:

(In practice, it's better to check a number of points along the line between B and E and make sure they're all land, because in some situations where B and E are straddling a long narrow bay the midpoint might be land with water to one side.)

Once a peninsula is detected, the next step is to determine where to anchor an island chain and which direction it should point.  This starts by placing a new point on the Middle point of the coastline between Beginning and End:

This point is the anchor for the island chain.  To measure how "pointy" the peninsula is, I create a line from the Middle point to the line between Beginning and End:

The ratio of the length of this line to the length of the line between Beginning and End is a measure of the peninsula's "pointiness".  Pointier peninsulas are better for anchoring island chains.  I don't want to put an island chain off every peninsula, so I scan the coastline, rank candidate locations based upon their pointiness, and then use the top candidates.

Here's an example island with the top four peninsulas highlighted:

You might note that the blue arrow isn't always perpendicular to the red dotted line.  This happens because "M" is at the midpoint of the coast line between "B" and "E".  Since the coastline is wiggly, the midpoint is not always exactly in the "middle".  In the example at the lower right part of the map you can see how sliding the red dotted line along the coast sometimes finds somewhat odd "peninsulas".  However, in general the algorithm works pretty well.

The best candidate on this map turns out to be the lower middle peninsula, so that is where I'll add an island chain.  To size the mask for this chain, I'll make the base of the wedge the length of the dotted red line, and I'll make the length of wedge about three times the length of the blue line.  On the following map I've filled in the mask with solid land to illustrate it's size and shape:

Keeping the mask no wider than the peninsula helps it look like a continuation of the peninsula.  The length is more a matter of taste -- choose what looks good to you.

Before creating the island range, I perturb the mask to make the regular shape less obvious:

And now I can generate the islands:

Here are a couple more examples:

In some cases, the archipelago extends into an area where there's already land -- an island off the coast, or another peninsula.  Adding a check that the mask area of the archipelago does not contain any land would prevent this from happening.  I choose to let it happen -- at best it looks like a small chain of land that connects two other land masses, which I like.

That covers the basic process Dragons Abound uses to create small archipelagos.

Time for A Stretch!

I've posted recently about how to generate mountain chains along fault lines, using a standard technique for generating mountains from noise functions, and then masking those mountains onto the map.  This allows me to create a range that looks like it follows a fault line where two tectonic plates have come together to create an upwelling of mountains:

If you examine the picture above, you'll see that it looks most like a bunch of cone shaped mountains  randomly overlapping within the area of the mountain range.  That's because I'm using a noise function that is symmetrical in both the x and y dimensions.  But in reality, the mountains within a range caused by a fault line have a different sort of structure.  Look at mountain ranges in this vintage postcard of Los Angeles:

(Click on any image to see it larger)

As this illustration shows, the mountains actually form ridges that run parallel to the fault line.  The two plates pushing up against one another create multiple folds in the earth.  I'd like to generate mountains that look more like this than my approach. To do this, I'm going to take advantage of an (underused) feature of gradient noise.  

Recall that currently the code in my mountain range generator looks something like this:

 for(var i=0;i<world.loc.length;i++) {
 [...]
    var val = pn.ridgedMF(x*4, x*4, 0.8, 3, 0.90) * height;

Here I'm stepping through the world, and in every mountain location adding 3 octaves of ridged multi-fractal noise with 90% persistence.  (The 0.8 is just a dummy value for the z dimension.)  I'm multiplying the world location by 4 before generating noise so that the frequency of the peaks in the mountains looks about right.  If I used a number less than 4, the peaks would too gradual and far apart.

 And likewise if I used a number bigger than 4, the peaks come too frequently:

The "right" number for this multiplier is dependent upon your coordinate system and is largely a matter of taste.

Typically, the same number is used for both dimensions to get a symmetrical noise.  If the numbers differ, then the noise gets "stretched" in one dimension or the other.  For example, if I set the multiplier for just the x dimension to 16, the mountains look like this:

As you can see, the noise is now varying faster in one dimension (the x dimension) than in the other dimension (the y dimension).  This has the useful effect of creating "ridges" in the noise.  In this case, the ridges run along the y dimension.  The way this mountain range is situated, I want the ridges to run the other way.   To get that, I have to reverse the multipliers so that the noise is varying faster in the y dimension:

And now the ridges run along the mountain range -- much like they do in the Santa Monica Mountains around Los Angeles.  This looks a little too chaotic to me -- I think it looks better if I multiply the y dimension by 8 and reduce the x dimension factor to 2.

If I have a mountain range running north-south, I can just flip the multipliers to create ridges running north-south:

And for any angle in between, I can just adjust the multipliers accordingly.... Right?

Uh, no.  Unfortunately, this technique can only create ridges along the x and y dimensions.  That's easy enough to see from the first image where the multipliers were the same.  That didn't create ridges running at 45 degrees, it created symmetrical peaks:

So it's going to take more work to get ridges at any angle.

Essentially, the solution is to create the mountain range along the x axis and then rotate and transform it into the proper position.  A quick trip to StackOverflow provides the equations for rotating a point (px,py) around a point (ox, oy) by a angle theta:

rx = cos(theta) * (px-ox) - sin(theta) * (py-oy) + ox
ry = sin(theta) * (px-ox) + cos(theta) * (py-oy) + oy

Fortunately there are a few simplifications that will make this easier.  First, since the noise function is infinite, the mountain doesn't actually have to be created along the x axis, but can be created along any line parallel to the x axis.  So it's only necessary to rotate the range, not transform its location -- which means I can choose any point (ox,oy) to anchor the rotation.  If I choose to rotate around the origin (i.e., (ox, oy) == (0, 0)), the rotation equations are simplified to this:

rx = cos(theta) * px - sin(theta) * py
ry = sin(theta) * px + cos(theta) * py

Second, it isn't necessary to create the range and then rotate it -- it suffices to do a "virtual" rotation and generate the heightmap for a location as if it were the rotated location.  In other words, I set the height of a location based upon the noise at the rotated location:

    height[x,y] = noise(rx, ry)

In pseudo-code the basics look like this:

    // Angle to rotate the ridges

    // (x0, y0) and (x1, y1) are two points that define the fault line
    var theta = -Math.atan2((y1-y0), (x1-x0));
    for(var i=0;i<world.loc.length;i++) {
       var x = world.loc[i].coords[0];
       var y = world.loc[i].coords[1];
       // First rotate x and y so that the ridges will run along
       // the length of the mountain chain.
       var rx = Math.cos(theta)*x - Math.sin(theta)*y;
       var ry = Math.sin(theta)*x + Math.cos(theta)*y;

       // Calculate the mountain height at (x,y) using (rx,ry)
       var val = pn.ridgedMF(rx*2, ry*8, 0.8, 3, 0.90) * height; 

       // (Code to calculate the mask omitted)

       world.loc[i].h += mask*val;

With this implemented,  I can now put a mountain range at any angle and have the ridges running along the range:

(The criss-cross patterns on the mountains are just an artifact of the fairly simple 3D display.)

The only remaining issue is the perturbation of the tent mask that creates the mountain range.  Recall that I perturb the tent mask with noise to give the mountain range a more organic shape instead of a ruler-straight rectangle.   To have the mountain ridges follow the same perturbed outline, I have to apply both the rotation transformation and the perturbation before setting the height.  Then I get ridge lines that follow the curves of the mountain range:

And now I have mountain ridges running organically along the underlying fault line.  This trick to create structure in the noise is pretty one-dimensional (so to speak) but it happens to be very handy for this one case of using noise to create mountain ranges.

All Good Things Must Come to An End

In my previous posting about mountain ranges along fault lines, I simplified things by assuming that the fault line and the mountain range atop it ran from one edge of the map to the other edge.  So there was no need to end the mountain range on the map.  But of course, that's something we might want to do, so in this post I'll describe my approach to that problem.

Recall that we placed the mountain range on the map by using a sort of tent-shaped mask along the fault line:

This mask was created by measuring the distance from the fault line, and making the mask high close to the line and dropping to zero at some distance away from the line.  So the mountains show through fully where the mask is high and disappear where the mask is low.

To make the mountains end on the map, we need to modify the mask so that it stops where we want the mountains to end.  And in fact, the function I provided in the previous posting had an option to treat the mask as a line segment, so that past the end of the line segment the "distance to the line" becomes the distance to the end point of the line segment, and the mask gets a curved end like this:

So we can use this version of the mask to end the mountains at the end of the line segment.  Here's an example:

This works, but in practice it's too abrupt.  The round end of the mountain range is also a little odd looking.  A map view makes this more apparent:

A more natural-looking termination is for the mountain range to trail off with a long, gradually narrowing end:

To do this, we need to end our tent mask with a long wedge shape rather than the half-circle of the simple mask.  We can achieve this adding ends to our tent mask that look like this:

(And now our mask really looks like a pup tent!)  Now when we draw our mountains inside this mask, they'll diminish to a pointy end:

With the addition of some perturbation to avoid arrow-straight mountain chains, the result can be very convincingly realistic:

The technique of using masks to control noise-based terrain is very powerful, because it permits you to create natural-looking terrain while simultaneously having exact control of how that terrain is placed on the map.

So how do you calculate the tent mask?  This is more complicated than the simpler mask used in the earlier posting.  In the version of the tent mask below, I also permit the main part of the tent to vary in width and height from one end to the other.  Altogether, the mask code looks like this:

//
//  Tent Mask
//
//  (x,y) -- the point we are masking
//  (x0, y0) -- one end of the main part of the tent
//   w0, h0  -- the width and height of the x0 end of the tent
//  taper0   -- a boolean that indicates whether the tent should taper past x0
//  l0       -- the length of the taper
//  (x1, y1) -- the other end of the main part of the tent
//   w1, h1  -- the width and height of the x1 end of the tent
//  taper1   -- a boolean that indicates whether the tent should taper past x1
//  l1       -- the length of the taper
//  
//
function tentMask(x, y, x0, y0, w0, h0, x1, y1, w1, h1, taper0, l0, taper1, l1) {
    // The math breaks down if x0 is the same as x1; rather than add that special
    // case, I simply tweak x1 to avoid the problem.
    if (x0 == x1) {
       x1 += 0.0001;
    };
    // x0 should be less than x1.  If not, swap the two points.
    if (x0 > x1) {
       var xt = x1, yt = y1, ht = h1, wt = w1;
       x1 = x0; y1 = y0; h1 = h0; w1 = w0;
       x0 = xt; y0 = yt; h0 = ht; w0 = wt;
    };
    //  (xp, yp) is going to be the point on the line (x0,y0)->(x1,y1) where the
    //  the perpendicular from (x, y) intersects.  If "t" is between 0 and 1, then
    //  the perpendicular hits between (x0,y0) and (x1,y1).  If not, then the
    //  perpendicular is off one end or the other.
    var d = (x1-x0)*(x1-x0)+(y1-y0)*(y1-y0);
    var d01 = Math.sqrt(d);
    var t = ((x-x0)*(x1-x0)+(y-y0)*(y1-y0))/d;
    var xp = x0 + t*(x1-x0);
    var yp = y0 + t*(y1-y0);
    //  If we're in the main part of the tent, or we're off and end without a
    //  taper, then our mask is just proportional to the distance of (x,y) from
    //  the line.
    if ((t >= 0 && t <= 1) || (t < 0 && !taper0) || (t > 1 && !taper1)) {
       // Distance from (x1,y1) to the point where the perpendicular meets the line
       var dp1 = Math.sqrt((xp-x1)*(xp-x1)+(yp-y1)*(yp-y1));
       // Height at the line
       var hp = (dp1/d01)*(h0-h1) + h1;
       // Width at this point of the line
       var wp = ((x1-xp)/(x1-x0))*(w0-w1)+w1;
       // Distance from the (x,y) to the line
       var lp = Math.sqrt((xp-x)*(xp-x)+(yp-y)*(yp-y));
       // if lp > wp then the point is outside the width of the
       // trapezoid, and the mask is 0
       if (lp > wp) return 0;
       // This is the height of the mask at (x,y)
       return (1-(lp/wp))*hp;
    };
    // Less than x0.  In this case, we want to use the distance away
    // from x0 and the taper length l0 to figure out the taper.
    if (t < 0 && taper0) {
       var dp0 = Math.sqrt((xp-x0)*(xp-x0)+(yp-y0)*(yp-y0));
       // If dp0 > l0 then we're off the end of the taper
       if (dp0 >= l0 || h0 == 0) return 0;
       // Height goes linear from h0 to 0
       var hp = (1-(dp0/l0))*h0/2 + h0/2;
       // Width goes linear from w0 to 0
       var wp = (1-(dp0/l0))*w0;
       // What's the width of x,y?
       var lp = Math.sqrt((xp-x)*(xp-x)+(yp-y)*(yp-y));
       // If lp > wp we're outside the arrow head
       if (lp > wp) return 0;
       return (1-(lp/wp))*hp;
    };
    // More than x1 -- the reverse case.
    if (t > 1 && taper1) {
       var dp1 = Math.sqrt((xp-x1)*(xp-x1)+(yp-y1)*(yp-y1));
       // If dp1 > l1 then we're off the end of the taper
       if (dp1 >= l1 || h0 == 0) return 0;
       // Height goes linear from h0 to 0
       var hp = (1-(dp1/l1))*h1;
       // Width goes linear from w0 to 0
       var wp = (1-(dp1/l1))*w1;
       // What's the width of x,y?
       var lp = Math.sqrt((xp-x)*(xp-x)+(yp-y)*(yp-y));
       // If lp > wp we're outside the arrow head
       if (lp > wp) return 0;
       return (1-(lp/wp))*hp;
    };
    // Otherwise masked out.
    return 0;
};

(Let me know if you like seeing code on the blog or not.)  The ends of the tent are really special cases of the middle of the tent, so this can probably be written more compactly as a conjunction of three masks, but it may be more understandable this way. 

This tent mask uses linear interpolation (the sides of the tent are straight lines).  A more complex mask might use a different easing function to create a more natural profile.  For map drawing this difference would probably not be apparent, but if you were using mountains as part of a 3D environment, this might be worth pursuing.

An Excursion: Lots of Faults

After my previous posting about creating mountain ranges based on fault lines, over on Reddit Bergasms wondered what a map would look like with lots of small faults.  That's an interesting question, so I took a small excursion to try it out.

Bergasm asked about small faults, meaning ones that start and end on the map.  In my previous posting, I assumed that the faults stretched all the way across the map.  In my next posting I'll talk in detail about how to terminate faults on the map, but for the moment you'll just have to assume that can be done.

To start off with, I generated a bunch of faults that were about 15% of the map size and randomly placed them on the map.  I turned off perturbation, reduced the spread, and didn't generate any mountains so the faults would be a little more obvious:

You can see a scattering of pill-shaped faults.  This is 10 faults. 

The isolated faults are a little jarring, but you can see in the lower right where faults are overlapping that the terrain is starting to look a little interesting.  Here's a map with 25 faults:

This is making some interesting shapes.  We can turn on perturbation to break up the regular shapes of the faults:

I think that makes a pretty interesting base terrain, at least as interesting as what we get from a noise function.  (It may be somewhat less "realistic", but that might be a feature for fantasy maps.)  I've seen some discussion of terrain generation using fault lines, but those algorithms typically use thousands of fault lines to generate a more jagged landscape.

Mountains don't work as well with these short faults.  If every fault generates mountains the map is nothing but mountains.

So we have to reduce the occurrence rate of mountains rather drastically.  But even so the mountains don't look very good because they are so short that they look like blobs rather than ranges:

There's also the problem that mountains stacked on top of mountains (as will often happen with lots of overlapping faults) quickly get unrealistic.  The central mountains here are 2X the height of "regular" mountains:

Overall I think I prefer the look of a map with a couple of longer faults:

But using many small faults to generate a base terrain as an alternate to a noise base seems promising.

An Excursion: Field Work

Headed out today to climb Bear Church Rock in the Shenandoah National Park.

 Apple included for scale.

Over the course of 4.5 hours, I determined quite conclusively that:  (1) Mountains are high, and (2) Much harder to walk over than to generate procedurally.  Further analysis is ongoing, or will be as soon as the Ibuprofen kicks in.

It's Not My Fault!

In a previous posting, I talked about how to generate mountains in a heightmap by using noise, and also about how to place mountains on a map so that they "make sense".  By that I meant that mountains should (1) only be on part of the map, (2) have a natural shape on the map, and (3) should transition smoothly into the rest of the terrain.  The approach I suggested was to use a mask, based either on noise or the already existing terrain (putting the mountains where the terrain was highest).  That yielded fairly natural-looking terrain:

In this posting, I'm going to talk about another approach for placing mountains that is based on the notion of fault lines.

A fault line is a fracture in the rock of the Earth's surface.  This happens for various reasons, but one cause is when two tectonic plates push past each other.  This can cause the rock of one plate to be thrust upward as the other plate burrows beneath the first plate.  Tectonic action of this sort created mountain ranges such as the Rocky Mountains in North America.  The very reason that we have notions like mountain "ranges" is because this mechanism is so primary in creating natural mountains.  So how can we implement a (simple) fault line model?

To start with, we can model a fault as a straight line cutting across the map from one edge to the other, forming a shallow peak at the line:

Implementing this is fairly straightforward:  Randomly create a line that crosses the map and then add height to the land based upon how far away the land is from the line.  For this, it's useful to know how to calculated the distance from a point to a line.  Here's an implementation in Javascript:

//
//  Distance from a Point to a Line Segment
//
//  This calculates the distance from (x,y) to the line defined by the
//  two points (x0, y0) and (x1, y1).  If segment is true, the line
//  is treated as a line segment.
//
function distanceFromPointToSegment(x, y, x0, y0, x1, y1, segment) {
    var d = (x1-x0)*(x1-x0)+(y1-y0)*(y1-y0);
    var t = ((x-x0)*(x1-x0)+(y-y0)*(y1-y0))/d;
    // If t < 0 || t > 1 we're off the end of the segment, and
    // the distance is just the distance to the proper endpoint. 
    if (segment && t < 0) {
    return Math.sqrt((x-x0)*(x-x0)+(y-y0)*(y-y0));
    };
    if (segment && t > 1) {
    return Math.sqrt((x-x1)*(x-x1)+(y-y1)*(y-y1));
    };
    // Otherwise calculate the perpendicular point and then
    // the distance to that point.
    var xp = x0 + t*(x1-x0);
    var yp = y0 + t*(y1-y0);
    return Math.sqrt((x-xp)*(x-xp)+(y-yp)*(y-yp));
}; 

The inverse distance can then be used to add a slope to the map on either side of the fault line, as in this exaggerated example:

The routine above can also treat the line as a line segment.  If we use this to add slope, it creates a rounded end at each end point of the line segment:

(We'll see why this is useful in a later post.)

A map can also have multiple fault lines, as in this example:

This gets busy and "unnatural" looking pretty quickly, although on a larger scale map it would make sense to have more fault lines.

One problem with these fault lines is that they're unnaturally straight and regular, as is evident if we map one:

Even with some noise added to break up the outlines, it's very evident where the fault lies (so to speak).

We can fix this by using noise to perturb the fault line.  Some tweaking is necessary to get appropriate noise parameters, but the result can be a very organic shape:

So we now have a natural looking fault line, created by drawing a perturbed line across the map and adding a gentle slope to either side.  Now we'd like to add a mountain range along the fault line.

We already know how to add mountains to an area by using a mask that fades out at the edges of the area.  What we need is a mask for the area around the fault line.  But that's just a steeper and narrower version of the fault itself.  Here I've overlayed a mask on top of the fault:

Where the mask is high, the mountains will show through fully; where it is low they will fade away.

However, you may have noticed that while the fault is a natural curve, the mountains mask is a straight line.  We correct that by apply the same perturbation we applied to the fault line to the mountains mask:

As you can see, the mountains mask now follows the same curve as the fault.  However, the mask has also become distorted: wider in some places, narrower in others, and with a bulge in the middle.  This is the natural consequence of using noise to perturb the mask.  The amount of noise can be tweaked to get a result that is acceptable -- you just need to make sure you use the same perturbation to both the fault line and the mountains mask.  The final step then is to generate mountains using noise and add the mountains to the map using the mountains mask as a filter:

The result is a pretty realistic-looking mountain "range".   Here's what that looks like in a map view:

 
There are a lot of tweaks you can apply to get different results, such as peakier mountains:

A more gradual transition into the surrounding terrain:

Multiple fault lines:

And so on.

There are a couple of shortcomings with this approach.  In real fault lines, the slope is steeper on one side and lower on the other.  That could be modeled with a mask that is steeper on one side than the other.  I also haven't talked about how to end a mountain range on the map, rather than running it from edge to edge.  In a later posting I'll talk about that, and how this technique can also be applied to generating island chains.

Profiling and Optimizing

Ugh.

I  had written a fairly long and detailed posting on profiling and optimizing the map generation program, but somehow during editing it I managed to delete everything.  I closed without saving so that I could go back to the previous version, but apparently Blogger had auto-saved in the interim.  What follows will be a somewhat shorter version of the original article.  Unfortunately, I lost the original profiling results, so the numbers below are only approximations.

A few days ago a comment on Reddit asked about the time required for various parts of the map generation.  Generally speaking I believe in writing code as simply as possible, and only worry about efficiency and optimization when it becomes a significant problem.  I find that much of my code development is very ephemeral -- for example, I wrote three or four different versions of erosion code before I settled on the current solution.  It would waste a lot of time to optimize code that might get replaced.  But the comment at least made me interested enough to turn on the Chrome profiler and see where the program spends the bulk of its time.  I generated this map:

This generation took about 65 seconds on my desktop machine.  Here are the routines that used at least 10% of the total execution time:

(39%) pickNextTerritoryLoc:  Territories (such as Kunapwknek in the map above) are created by placing a capitol city at a good location and then growing the territory outward.  This routine is responsible for picking the next location.  It scans the map looking for locations adjacent to the territory, evaluates each location for desirability, and returns the best location it found.

(20%) drawLabels:  This is the routine responsible for labeling the map.  The map generator uses a method called "force layout" to place the labels.  Most of the time is spent in a function that finds an open spot for territory labels, and in a function that detects when bounding boxes collide as labels are moved around.

(20%) neighborsWithin:  This takes a location and finds all the neighbors within a given distance.  On a grid map this would be straightforward, but on a Voronoi grid it is more challenging because each location is irregularly sized and can have different numbers of neighbors.

(15%) windModel:  The wind model is responsible for propagating trade winds across the map.  Because winds can be affected by the rising and falling land, and because the Voronoi grid is irregular, this is a complex iterative process.

For different reasons, I'm not interested in trying to optimize drawLabels or windModel right now.  I may replace drawLabels at some point with a simulated annealing algorithm, so it may not make much sense to sink a lot of time into optimizing it.  And I've already spent some time making the windModel faster, so there's probably not a lot of easier improvement available there.  But pickNextTerritoryLoc and neighborsWithin are worth examining to see if there are any easy gains to be had.

Examining the code for pickNextTerritoryLoc reveals two possible avenues for optimizing.  The first opportunity is the function for evaluating locations.  It is called every time a location is considered, but most of the evaluation is static and could be computed once before the territories are grown.  The second opportunity is to improve how the routine identifies candidate locations.  Currently the map is scanned to identify all the territory's neighbors.  This is inefficient when the routine is called many times.

(As a side note, when you're implementing optimizations, it's best to work in small steps and check after every step that the program behavior hasn't changed.  This guards against introducing errors and gives you a chance to monitor the function efficiency.  Sometimes an early change saves so much time that doing the rest is unnecessary.  In this case, I began by splitting the evaluation function into static and dynamic components.  In doing so, I discovered a bug in the original function that I had to correct before attempting the optimization.)

Splitting the evaluation function into static and dynamic pieces, and only calling the static portion once yielded about a 2% gain.  As it happens, the most time-consuming part of the evaluation had to remain in the dynamic piece, and that limited the improvement.  The second optimization was to have each territory maintain a list of neighboring locations and incrementally update that list as the territory grew, rather than repeatedly scanning the map.  This was much more fruitful, and dropped the percentage of total time taken by this function from 37% to 5%.  Good enough!

Moving on to neighborsWithin, an examination of the code didn't reveal any obvious optimizations.  I attempted an optimization that involved switching from arrays of neighbors to a hash table, but partway through rewriting the code realized that this approach wouldn't work and backed out the changes.  I also thought about using a quadTree representation to speed up the search.  This would probably work (and is on my To-Do list, because it would be potentially useful in many places) but it would be more work than I was willing to spend at the moment.

In the end, I fell back on the old optimization trick of trading off space for time.  The real problem with neighborsWithin is that it gets called many times for each location on the map.  A simple optimization strategy in this situation is to cache the function result.  The first time the function is asked for the neighbors within some distance of a location it calculates the result, but it also saves the result.  The next time the function is asked the same question, it uses the stored result rather than calculate it again.  (Obviously this only works if the answer doesn't change, or if you can tell whether it has changed and needs re-calculation.)  This change dropped neighborsWithin from 20% to ~4.5% of the computing time.  Again, that was a nice gain for very little effort!

These two optimizations dropped the overall map generation time from 65 seconds to 44 seconds.  After all that, here's the approximate times spent in different phases of the map generation and display:

(27%) Generating the World (physical geography)
  (17%) Wind Model
  (4%) Add a Central Island
  (4%) Creating Biomes
(4%) Generating the Culture (countries, etc.)
(57%) Drawing the Map
  (55%) Drawing the Labels
     (40%) Finding an Open Spot for Territory Labels
     (16%) Executing the Force Layout

The glaring problem here is that 40% of the total execution time is now spent on the relatively minor task of finding an open spot for the territory labels!  As I mentioned above, I'm reluctant to spend too much time on optimizing the label layout because I'm already half-convinced that I should switch to a simulated annealing approach, but this may be worth at least a preliminary look.