From 9b9a471aa80c516cb06785f6246b74ec2f0434f6 Mon Sep 17 00:00:00 2001
From: Starbeamrainbowlabs <sbrl@starbeamrainbowlabs.com>
Date: Tue, 9 Jun 2020 01:21:32 +0100
Subject: [PATCH] Implement initial (untested) convolution system. Next we need
 to implement a worldedit function to handle fetching the manip data,
 calculating the heightmap, pushing it through this convolutional system, and
 saving the changes back again.

---
 .../lib/convolution/convolution.lua           |  6 ++
 .../lib/convolution/convolve.lua              | 43 +++++++++
 .../lib/convolution/kernel_gaussian.lua       | 51 ++++++++++
 .../lib/convolution/kernels.lua               | 94 +++++++++++++++++++
 worldeditadditions/utils/numbers.lua          |  6 ++
 worldeditadditions/utils/strings.lua          | 26 ++++-
 6 files changed, 224 insertions(+), 2 deletions(-)
 create mode 100644 worldeditadditions/lib/convolution/convolution.lua
 create mode 100644 worldeditadditions/lib/convolution/convolve.lua
 create mode 100644 worldeditadditions/lib/convolution/kernel_gaussian.lua
 create mode 100644 worldeditadditions/lib/convolution/kernels.lua

diff --git a/worldeditadditions/lib/convolution/convolution.lua b/worldeditadditions/lib/convolution/convolution.lua
new file mode 100644
index 0000000..b68663f
--- /dev/null
+++ b/worldeditadditions/lib/convolution/convolution.lua
@@ -0,0 +1,6 @@
+worldeditadditions.conv = {}
+
+dofile(worldeditadditions.modpath.."/lib/conv/kernels.lua")
+dofile(worldeditadditions.modpath.."/lib/conv/kernel_gaussian.lua")
+
+dofile(worldeditadditions.modpath.."/lib/conv/convolve.lua")
diff --git a/worldeditadditions/lib/convolution/convolve.lua b/worldeditadditions/lib/convolution/convolve.lua
new file mode 100644
index 0000000..c5d1c71
--- /dev/null
+++ b/worldeditadditions/lib/convolution/convolve.lua
@@ -0,0 +1,43 @@
+
+--[[
+Convolves over a given 2D heightmap with a given matrix.
+Note that this *mutates* the given heightmap.
+Note also that the dimensions of the matrix must *only* be odd.
+@param	{number[]}			heightmap		The 2D heightmap to convolve over.
+@param	{[number,number]}	heightmap_size	The size of the heightmap as [ height, width ]
+@param	{number[]}			matrix			The matrix to convolve with.
+@param	{[number, number]}	matrix_size		The size of the convolution matrix as [ height, width ]
+]]--
+function worldeditadditions.conv.convole(heightmap, heightmap_size, matrix, matrix_size)
+	if matrix_size[0] % 2 ~= 1 or matrix_size[1] % 2 ~= 1 then
+		return false, "Error: The matrix size must contain only odd numbers (even number detected)"
+	end
+	
+	local border_size = {
+		(matrix_size[0]-1) / 2,	-- height
+		(matrix_size[1]-1) / 2	-- width
+	}
+	
+	-- Convolve over only the bit that allows us to use the full convolution matrix
+	for y = heightmap_size[0] - border_size[0], border_size[0], -1 do
+		for x = heightmap_size[1] - border_size[1], border_size[1], -1 do
+			local total = 0
+			
+			for my = matrix_size[0], 0, -1 do
+				for mx = matrix_size[1], 0, -1 do
+					local mi = (my * matrix_size[1]) + mx
+					local cy = y + (my - border_size[0])
+					local cx = x + (mx - border_size[1])
+					
+					local i = (cy * heightmap_size[1]) + cx
+					
+					total = total + matrix[mi] * heightmap[i]
+				end
+			end
+			
+			heightmap[(y * heightmap_size[1]) + x] = total
+		end
+	end
+	
+	return true, heightmap
+end
diff --git a/worldeditadditions/lib/convolution/kernel_gaussian.lua b/worldeditadditions/lib/convolution/kernel_gaussian.lua
new file mode 100644
index 0000000..d4025a9
--- /dev/null
+++ b/worldeditadditions/lib/convolution/kernel_gaussian.lua
@@ -0,0 +1,51 @@
+-- Ported from Javascript by Starbeamrainbowlabs
+-- Original source: https://github.com/sidorares/gaussian-convolution-kernel/
+-- From 
+-- the code is taken from https://github.com/mattlockyer/iat455/blob/6493c882f1956703133c1bffa1d7ee9a83741cbe/assignment1/assignment/effects/blur-effect-dyn.js
+-- (c) Matt Lockyer, https://github.com/mattlockyer
+
+-- hypotenuse has moved to utils/numbers.lua
+
+--[[
+ * Generates a kernel used for the gaussian blur effect.
+ *
+ * @param dimension is an odd integer
+ * @param sigma is the standard deviation used for our gaussian function.
+ *
+ * @returns an array with dimension^2 number of numbers, all less than or equal
+ *	 to 1. Represents our gaussian blur kernel.
+]]--
+function worldeditadditions.conv.kernel_gaussian(dimension, sigma)
+	if not (dimension % 2) or math.floor(dimension) ~= dimension or dimension < 3 then
+		return false, "The dimension must be an odd integer greater than or equal to 3"
+	end
+	local kernel = {};
+
+	local two_sigma_square = 2 * sigma * sigma;
+	local centre = (dimension - 1) / 2;
+	
+	local sum = 0
+	for i = 0, dimension-1 do
+		for j = 0, dimension-1 do
+			local distance = worldeditadditions.hypotenuse(i, j, centre, centre)
+
+			-- The following is an algorithm that came from the gaussian blur
+			-- wikipedia page [1].
+			--
+			-- http://en.wikipedia.org/w/index.php?title=Gaussian_blur&oldid=608793634#Mechanics
+			local gaussian = (1 / math.sqrt(
+				math.pi * two_sigma_square
+			)) * math.exp((-1) * (math.pow(distance, 2) / two_sigma_square));
+			
+			sum = sum + gaussian
+			kernel[i*dimension + j] = gaussian
+		end
+	end
+	
+	-- Returns the unit vector of the kernel array.
+	for k,v in pairs(kernel) do
+		kernel[k] = kernel[k] / sum
+	end
+	
+	return kernel
+end
diff --git a/worldeditadditions/lib/convolution/kernels.lua b/worldeditadditions/lib/convolution/kernels.lua
new file mode 100644
index 0000000..3f89144
--- /dev/null
+++ b/worldeditadditions/lib/convolution/kernels.lua
@@ -0,0 +1,94 @@
+--- Creates a (normalised) box convolutional kernel.
+-- Larger box kernels will obviously be slower, but will produce a more blurred
+-- effect (i.e. smoother terrain).
+-- @param	width	number	The width of the kernel.
+-- @param	height	number	The height of the kernel.
+-- @return	The resulting kernel as a ZERO-indexed list of numbers.
+function worldeditadditions.conv.kernel_box(width, height)
+	local result = {}
+	local total = 0
+	for y = 0, height do
+		for x = 0, width do
+			result[(y*width) + x] = 1
+			total = total + 1
+		end
+	end
+	-- Ensure that everything sums up to 1
+	for i = 0, #result do
+		result[i] = result[i] / total
+	end
+	return result
+end
+
+--- Computes the Lth line of Pascal's triangle.
+-- More information: https://en.wikipedia.org/wiki/Pascal%27s_triangle
+-- There are probably more efficient ways to it that don't repeat themselves as
+-- much, but this is my solution.
+-- @param l				number	The 1-indexed row of Pascal's Triangle to return.
+-- @return number[]		A ZERO-indexed list of numbers in the specified row of Pascal's Triangle.
+local function pascal(l)
+	local prev = {}
+	prev[0] = 1
+	if l == 1 then return prev end
+	prev[1] = 1
+	if l == 2 then return prev end
+	local length_last = 2
+	
+	for n=3, l do
+		local next = {}
+		for i=0, length_last do
+			if i == 0 or i == length_last then
+				next[i] = 1
+			else
+				next[i] = prev[i - 1] + prev[i]
+			end
+		end
+		prev = next
+		length_last = length_last + 1
+	end
+	return prev
+end
+
+--- Creates a pascal convolutional kernel.
+-- Larger pascal kernels will obviously be slower, but will produce a more blurred
+-- effect (i.e. smoother terrain).
+-- @param	width			number	The width of the kernel.
+-- @param	height			number	The height of the kernel.
+-- @param	normalise=true	bool	Whether to normalise the resulting kernel (default: true)
+-- @return	The resulting kernel as a ZERO-indexed list of numbers.
+function worldeditadditions.conv.kernel_pascal(width, height, normalise)
+	if normalise == nil then normalise = true end
+	
+	local result = {}
+	local pascal_width = width
+	local height_middle = ((height - 2) / 2)
+	local total = 0
+	for y = 0, height-1 do
+		local pascal_numbers = pascal(pascal_width)
+		local pascal_start = (pascal_width - width) / 2
+		for x = 0, width - 1 do
+			result[(y*width) + x] = pascal_numbers[pascal_start + x]
+			total = total + pascal_numbers[pascal_start + x]
+		end
+		
+		if y > height_middle then pascal_width = pascal_width - 2
+		else pascal_width = pascal_width + 2 end
+	end
+	
+	if normalise then
+		for k,v in pairs(result) do
+			result[k] = result[k] / total
+		end
+	end
+	
+	return result
+end
+
+-- print("box, 5x5")
+-- print_2d(kernel_box(5, 5), 5)
+-- print("pascal, 5x5")
+-- print_2d(kernel_pascal(5, 5), 5)
+-- print("pascal, 7x7")
+-- print_2d(kernel_pascal(7, 7), 7)
+-- print("pascal, 9x9")
+-- print_2d(kernel_pascal(9, 9), 9)
diff --git a/worldeditadditions/utils/numbers.lua b/worldeditadditions/utils/numbers.lua
index c54bf54..06a9d63 100644
--- a/worldeditadditions/utils/numbers.lua
+++ b/worldeditadditions/utils/numbers.lua
@@ -3,3 +3,9 @@ function worldeditadditions.round(num, numDecimalPlaces)
   local mult = 10^(numDecimalPlaces or 0)
   return math.floor(num * mult + 0.5) / mult
 end
+
+function worldeditadditions.hypotenuse(x1, y1, x2, y2)
+	local xSquare = math.pow(x1 - x2, 2);
+	local ySquare = math.pow(y1 - y2, 2);
+	return math.sqrt(xSquare + ySquare);
+end
diff --git a/worldeditadditions/utils/strings.lua b/worldeditadditions/utils/strings.lua
index 844f45a..5086624 100644
--- a/worldeditadditions/utils/strings.lua
+++ b/worldeditadditions/utils/strings.lua
@@ -73,10 +73,31 @@ function worldeditadditions.str_padstart(str, len, char)
     return string.rep(char, len - #str) .. str
 end
 
+--- Equivalent to string.startsWith in JS
+-- @param	str		string	The string to operate on
+-- @param	start	number	The start string to look for
+-- @returns	bool	Whether start is present at the beginning of str
 function worldeditadditions.string_starts(str,start)
    return string.sub(str,1,string.len(start))==start
 end
 
+--- Prints a 2d array of numbers formatted like a JS TypedArray (e.g. like a manip node list or a convolutional kernel)
+-- In other words, the numbers should be formatted as a single flat array.
+-- @param	tbl		number[]	The ZERO-indexed list of numbers
+-- @param	width	number		The width of 2D array.
+function worldeditadditions.print_2d(tbl, width)
+	local next = {}
+	for i=0, #tbl do
+		table.insert(next, tbl[i])
+		if #next == width then
+			print(table.concat(next, "\t"))
+			next = {}
+		end
+	end
+end
+
+
+
 --- Turns an associative node_id → count table into a human-readable list.
 -- Works well with worldeditadditions.make_ascii_table().
 function worldeditadditions.node_distribution_to_list(distribution, nodes_total)
@@ -127,8 +148,9 @@ function worldeditadditions.make_ascii_table(data, total)
 	return table.concat(result, "\n")
 end
 
-
-
+--- Parses a list of strings as a list of weighted nodes - e.g. like in the //mix command.
+-- @param	parts	string[]	The list of strings to parse (try worldeditadditions.split)
+-- @returns	table	A table in the form node_name => weight.
 function worldeditadditions.parse_weighted_nodes(parts)
 	local MODE_EITHER = 1
 	local MODE_NODE = 2