# HG changeset patch
# User Teemu Piippo <teemu.s.piippo@gmail.com>
# Date 1656853162 -10800
# Node ID 180072db4a834184d25d1a78efc7dc114220302c
# Parent  af6633412a6cd72dfef6df6879f0e0d3b60e02da
Fix polygonize not finding the normal properly

diff -r af6633412a6c -r 180072db4a83 src/geometry.h
--- a/src/geometry.h	Sun Jul 03 14:35:06 2022 +0300
+++ b/src/geometry.h	Sun Jul 03 15:59:22 2022 +0300
@@ -198,3 +198,26 @@
 
 glm::vec3 pointOnCurve(const BezierCurve& curve, float t);
 glm::vec3 derivativeOnCurve(const BezierCurve& curve, float t);
+
+template<typename Iter>
+std::optional<glm::mat3> calculateNormal(Iter&& begin, Iter&& end)
+{
+	const long int n = end - begin;
+	std::optional<glm::mat3> result;
+	for (long int i = 0; i < n; ++i) {
+		const glm::vec3& v1 = *(begin + (i + n - 1) % n);
+		const glm::vec3& v2 = *(begin + i);
+		const glm::vec3& v3 = *(begin + (i + 1) % n);
+		const glm::vec3 xvec = v1 - v2;
+		const glm::vec3 yvec = v3 - v2;
+		constexpr float threshold = 1e-6f;
+		if (glm::length(xvec) > threshold and glm::length(yvec) > threshold) {
+			const glm::vec3 zvec = glm::cross(glm::normalize(xvec), glm::normalize(yvec));
+			if (glm::length(zvec) > threshold) {
+				result = {xvec, yvec, zvec};
+				break;
+			}
+		}
+	}
+	return result;
+}
diff -r af6633412a6c -r 180072db4a83 src/triangulate.cpp
--- a/src/triangulate.cpp	Sun Jul 03 14:35:06 2022 +0300
+++ b/src/triangulate.cpp	Sun Jul 03 15:59:22 2022 +0300
@@ -134,28 +134,29 @@
 	std::vector<glm::vec3>::const_iterator end)
 {
 	std::vector<PlainPolygonElement> result;
-	const glm::vec3 xvec = *begin - *(begin + 1);
-	const glm::vec3 yvec = *(begin + 2) - *(begin + 1);
-	const glm::vec3 zvec = glm::cross(xvec, yvec);
-	const glm::mat3 inverseMatrix = glm::inverse(glm::mat3{xvec, yvec, zvec});
-	std::vector<std::vector<glm::vec3>> polygons{1};
-	std::vector<glm::vec3>& polygon2d = polygons.back();
-	polygon2d.reserve(static_cast<std::size_t>(end - begin));
-	for (std::vector<glm::vec3>::const_iterator it = begin; it != end; ++it) {
-		polygon2d.push_back(inverseMatrix * glm::vec4{*it, 1});
-	}
-	const std::vector<triangulate::index_t> indices = mapbox::earcut<triangulate::index_t>(polygons);
-	triangulate::MergedTriangles mergedTriangles = triangulate::mergeTriangles(indices, begin);
-	for (Quadrilateral& quad : mergedTriangles.quadrilaterals) {
-		result.push_back(quad);
-	}
-	for (std::size_t i = 0; i < indices.size(); i += 3) {
-		if (not mergedTriangles.cutTriangles.contains({i})) {
-			Triangle tri = triangle(
+	// Transform the polygon into XY plane
+	opt<glm::mat3> normal = calculateNormal(begin, end);
+	if (normal.has_value()) {
+		const glm::mat3 normalInverse = glm::inverse(*normal);
+		std::vector<std::vector<glm::vec3>> polygons{1};
+		std::vector<glm::vec3>& polygon2d = polygons.back();
+		polygon2d.reserve(static_cast<std::size_t>(end - begin));
+		for (std::vector<glm::vec3>::const_iterator it = begin; it != end; ++it) {
+			polygon2d.push_back(normalInverse * glm::vec4{*it, 1});
+		}
+		const std::vector<triangulate::index_t> indices = mapbox::earcut<triangulate::index_t>(polygons);
+		triangulate::MergedTriangles mergedTriangles = triangulate::mergeTriangles(indices, begin);
+		for (Quadrilateral& quad : mergedTriangles.quadrilaterals) {
+			result.push_back(quad);
+		}
+		for (std::size_t i = 0; i < indices.size(); i += 3) {
+			if (not mergedTriangles.cutTriangles.contains({i})) {
+				Triangle tri = triangle(
 				*(begin + indices[i]),
 				*(begin + indices[i + 1]),
 				*(begin + indices[i + 2]));
-			result.push_back(tri);
+				result.push_back(tri);
+			}
 		}
 	}
 	return result;