/**
 * Usage:
 *   towns = GSTownList();
 *   towns.Valuate(GSTown.GetLocation);
 *   neighbours = Delaunay(towns);
 */
require("graph.nut");
require("algebra.nut");

function CircleCondition(graph, edge, triangles) { // returns true, if the circumcircle of one triangle contains a point of another triangle
	// calculate perpendicular bisectors
	local perp1 = Perpendicular([[edge[0].x, edge[0].y], [edge[1].x, edge[1].y]]);
	local perp2 = Perpendicular([[edge[0].x, edge[0].y], [triangles[0].x, triangles[0].y]]);
	local center = Intersection2(perp1, perp2);
	if (center == null) return true; // should never happen, but apparently it does
	return Dist2sq(center, [edge[0].x, edge[0].y]) > Dist2sq(center, [triangles[1].x, triangles[1].y]);
}

function Flip(graph, edge) {
	local triangles = graph.FindTriangles(edge[0], edge[1]);
	if (triangles.len() < 2) return;
	if (!CircleCondition(graph, edge, triangles)) return;
	//GSLog.Info("Flip: <" + edge[0] + ", " + edge[1] + "> - <" + triangles[0] + ", " + triangles[1] + ">");
	graph.RemoveEdge(edge[0], edge[1]);
	graph.AddEdge(triangles[0], triangles[1]);
	graph.flips++;
	// recurse
	if (graph.HasEdge(edge[0], triangles[0])) Flip(graph, [edge[0], triangles[0]]);
	if (graph.HasEdge(edge[0], triangles[1])) Flip(graph, [edge[0], triangles[1]]);
	if (graph.HasEdge(edge[1], triangles[0])) Flip(graph, [edge[1], triangles[0]]);
	if (graph.HasEdge(edge[1], triangles[1])) Flip(graph, [edge[1], triangles[1]]);
}

function Delaunay(town_list) { // a list valuated by tileindex
	GSLog.Info("Delaunay Triangulation starting");
	local neighbours = Graph(); // undirected graph
	local projections = [];     // directed edges
	local a = neighbours.AddNode(town_list.Begin()); // initialize with edge between first two towns
	//GSLog.Info("Town: " + a + " " + GSTown.GetName(GSTile.GetClosestTown(a.ti)));
	if (town_list.IsEnd()) return neighbours; // edge case: one town
	local b = neighbours.AddNode(town_list.Next());
	//GSLog.Info("Town: " + b + " " + GSTown.GetName(GSTile.GetClosestTown(b.ti)));
	neighbours.AddEdge(a, b);
	projections.append([a,b]);
	if (a.y != b.y) projections.append([b,a]);
	town_list.RemoveValue(a.ti);
	town_list.RemoveValue(b.ti);
	foreach (g, _ in town_list) {
		local t = neighbours.AddNode(g);
		//GSLog.Info("Town: " + t + " " + GSTown.GetName(GSTile.GetClosestTown(t.ti)));
		local plane = t.y;
		local new_projections = [];
		local last = null;
		local parallel_node = null;
		// TODO: binary search of matching projections
		foreach (edge in projections) {
			local point = Intersection([[edge[0].x, edge[0].y], [edge[1].x, edge[1].y]], [[0,plane],[GSMap.GetMapSizeX(),plane]])
			/*if (point != null) {
				GSLog.Info("Intersection: <" + edge[0] + ", " + edge[1] + "> - (" + point[0] + ", " + point[1] + ") " + point[2]);
			} else {
				GSLog.Info("Intersection: <" + edge[0] + ", " + edge[1] + "> - parallel");
			}*/
			if (
				(point == null) && (edge[0].y != plane) || // parallel
				(point != null) && (
					(point[2] <= 0) && (t.x <= point[0]) ||
					(point[2] >= 0) && (t.x >= point[0])
				)
			) {
				parallel_node = null;
				neighbours.AddEdge(edge[0], t);
				//GSLog.Info("Add: <" + edge[0] + ", " + t + ">");
				neighbours.AddEdge(edge[1], t);
				//GSLog.Info("Add: <" + edge[1] + ", " + t + ">");
				if (neighbours.HasEdge(edge[0], edge[1])) Flip(neighbours, edge);
				if (neighbours.HasEdge(edge[0], t))       Flip(neighbours, [edge[0], t]);
				if (last == null) {
					//GSLog.Info("Insert: <" + edge[0] + ", " + t + ">");
					new_projections.append([edge[0], t]);
				}
				last = [t, edge[1]];
			} else if (point == null) /*&& (edge[0].y == plane)*/ {
				parallel_node = edge[1];
				new_projections.append(edge);
			} else {
				parallel_node = null;
				if (last != null) {
					//GSLog.Info("Insert: <" + last[0] + ", " + last[1] + ">");
					new_projections.append(last);
					last = null;
				}
				if ((point[0] > 0) && (point[0] < GSMap.GetMapSizeX())) {
					//GSLog.Info("Keep: <" + edge[0] + ", " + edge[1] + ">");
					new_projections.append(edge);
				}/* else {
					GSLog.Info("Skip: <" + edge[0] + ", " + edge[1] + ">");
				}*/
			}
		}
		if (parallel_node != null) {
			neighbours.AddEdge(parallel_node, t);
			//GSLog.Info("Add: <" + parallel_node + ", " + t + ">");
			last = [parallel_node, t];
		}
		if (last != null) {
			//GSLog.Info("Insert: <" + last[0] + ", " + last[1] + ">");
			new_projections.append(last);
		}
		projections = new_projections;
		/*//local s="";
		local prev = null;
		foreach (i, edge in projections) {
			if ((prev != null) && (prev[1] != edge[0])) GSLog.Error("Error: inconsistent projections");
			prev = edge;
			//if (i>0) s+=", ";
			//s += edge[0];
		}*/
		//GSLog.Info("Graph: " + neighbours);
		//GSLog.Info("Projections: ["+s+"]");
	}
	GSLog.Info("Delaunay Triangulation finished (" + neighbours.nodes.len() + " Towns, " + neighbours.flips + " Flips)");
	return neighbours;
}