commit 0aadd793f5f2ec9c3460428e57fd5fd20eadd14a
Author: Marc Lepage <mlepage@antimeta.com>
Date:   Sat Oct 22 18:06:29 2011 -0400

    Initial commit.

diff --git a/COPYRIGHT b/COPYRIGHT
new file mode 100644
index 0000000..c165c29
--- /dev/null
+++ b/COPYRIGHT
@@ -0,0 +1,19 @@
+Copyright (C) 2011 Marc Lepage
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in
+all copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+THE SOFTWARE.
diff --git a/README b/README
new file mode 100644
index 0000000..8d8ebfb
--- /dev/null
+++ b/README
@@ -0,0 +1,66 @@
+Heightmap module
+by Marc Lepage
+
+
+OVERVIEW
+
+The heightmap module uses the diamond-square algorithm to generate cloud or
+plasma fractal heightmaps which can be used for terrain.
+
+
+USAGE
+
+-- import module
+require "heightmap"
+
+-- create 256x256 heightmap
+map = heightmap.create(256, 256)
+
+-- examine each height value
+for x = 0, map.w do
+    for y = 0, map.y do
+        print(map[x][y])
+    end
+end
+
+-- define a custom height function (reusing the default but scaling it)
+function f(map, x, y, d, h)
+    return 2 * heightmap.defaultf(map, x, y, d, h)
+end
+
+-- use it to create a new heightmap
+map = heightmap.create(256, 256, f)
+
+
+HOW IT WORKS
+
+The heightmap must be a square the size of a power of two, plus one, so that
+it can be evenly divided. For example, 4x4 cells will require 5x5 vertices.
+
+First the four corners are seeded with a random value (C).
+
+Then each square is used to set the value of its center (S) based on the
+average of its four corners (plus some randomness).
+
+Then each diamond is used to set the value of its center (D) based on the
+average of its four points (plus some randomness).
+
+The square and diamond steps continue until all values have been set:
+
+      4     S 2    D 2    S 1    D 1
+    C...C  c...c  c.D.c  c.d.c  cDdDc
+    .....  .....  .....  .S.S.  DsDsD
+    .....  ..S..  D.s.D  d.S.d  dDsDd
+    .....  .....  .....  .S.S.  DsDsD
+    C...C  c...c  c.D.c  c.d.c  cDdDc
+
+The default height function randomly displaces values by up to +/- 0.5 of the
+step size. So above, the corners will be from -2 to +2, the center will be
+the mean of the corners randomly displaced from -1 to +1, and so on.
+
+
+RESOURCES
+
+http://en.wikipedia.org/wiki/Diamond-square_algorithm
+http://en.wikipedia.org/wiki/Heightmap
+http://en.wikipedia.org/wiki/Fractal_landscape
diff --git a/heightmap.lua b/heightmap.lua
new file mode 100644
index 0000000..e132283
--- /dev/null
+++ b/heightmap.lua
@@ -0,0 +1,106 @@
+-- Heightmap module
+-- Copyright (C) 2011 Marc Lepage
+
+local max, random = math.max, math.random
+
+module(...)
+
+-- Find power of two sufficient for size
+local function pot(size)
+	local pot = 2
+	while true do
+		if size <= pot then return pot end
+		pot = 2*pot
+	end
+end
+
+-- Create a table with 0 to n zero values
+local function tcreate(n)
+	local t = {}
+	for i = 0, n do t[i] = 0 end
+	return t
+end
+
+-- Square step
+-- Sets map[x][y] from square of radius d using height function f
+local function square(map, x, y, d, f)
+	local sum, num = 0, 0
+	if 0 <= x-d then
+		if   0 <= y-d   then sum, num = sum + map[x-d][y-d], num + 1 end
+		if y+d <= map.h then sum, num = sum + map[x-d][y+d], num + 1 end
+	end
+	if x+d <= map.w then
+		if   0 <= y-d   then sum, num = sum + map[x+d][y-d], num + 1 end
+		if y+d <= map.h then sum, num = sum + map[x+d][y+d], num + 1 end
+	end
+	map[x][y] = f(map, x, y, d, sum/num)
+end
+
+-- Diamond step
+-- Sets map[x][y] from diamond of radius d using height function f
+local function diamond(map, x, y, d, f)
+	local sum, num = 0, 0
+	if   0 <= x-d   then sum, num = sum + map[x-d][y], num + 1 end
+	if x+d <= map.w then sum, num = sum + map[x+d][y], num + 1 end
+	if   0 <= y-d   then sum, num = sum + map[x][y-d], num + 1 end
+	if y+d <= map.h then sum, num = sum + map[x][y+d], num + 1 end
+	map[x][y] = f(map, x, y, d, sum/num)
+end
+
+-- Diamond square algorithm generates cloud/plasma fractal heightmap
+-- http://en.wikipedia.org/wiki/Diamond-square_algorithm
+-- Size must be power of two
+-- Height function f must look like f(map, x, y, d, h) and return h'
+local function diamondsquare(size, f)
+	-- create map
+	local map = { w = size, h = size }
+	for c = 0, size do map[c] = tcreate(size) end
+	-- seed four corners
+	local d = size
+	map[0][0] = f(map, 0, 0, d, 0)
+	map[0][d] = f(map, 0, d, d, 0)
+	map[d][0] = f(map, d, 0, d, 0)
+	map[d][d] = f(map, d, d, d, 0)
+	d = d/2
+	-- perform square and diamond steps
+	while 1 <= d do
+		for x = d, map.w-1, 2*d do
+			for y = d, map.h-1, 2*d do
+				square(map, x, y, d, f)
+			end
+		end
+		for x = d, map.w-1, 2*d do
+			for y = 0, map.h, 2*d do
+				diamond(map, x, y, d, f)
+			end
+		end
+		for x = 0, map.w, 2*d do
+			for y = d, map.h-1, 2*d do
+				diamond(map, x, y, d, f)
+			end
+		end
+		d = d/2
+	end
+	return map
+end
+
+-- Default height function
+-- d is depth (from size to 1 by powers of two)
+-- h is mean height at map[x][y] (from square/diamond of radius d)
+-- returns h' which is used to set map[x][y]
+function defaultf(map, x, y, d, h)
+	return h + (random()-0.5)*d
+end
+
+-- Create a heightmap using the specified height function (or default)
+-- map[x][y] where x from 0 to map.w and y from 0 to map.h
+function create(width, height, f)
+	f = f and f or defaultf
+	-- make heightmap
+	local map = diamondsquare(pot(max(width, height)), f)
+	-- clip heightmap to desired size
+	for x = 0, map.w do for y = height+1, map.h do map[x][y] = nil end end
+	for x = width+1, map.w do map[x] = nil end
+	map.w, map.h = width, height
+	return map
+end