smallest circumcircle and circumcenter
diff --git a/CMakeLists.txt b/CMakeLists.txt
index c8e8e13..c12ca18 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -1,6 +1,6 @@
cmake_minimum_required(VERSION 3.0.0)
project(delaunator VERSION 0.1.0)
-set (CMAKE_CXX_STANDARD 11)
+set (CMAKE_CXX_STANDARD 17)
if (NOT EXISTS "${PROJECT_SOURCE_DIR}/includes")
execute_process(COMMAND bash "-c" "(cd ${PROJECT_SOURCE_DIR} && ./fetch-includes.sh)")
endif()
diff --git a/src/delaunator.cpp b/src/delaunator.cpp
index e3d37fe..c05b127 100644
--- a/src/delaunator.cpp
+++ b/src/delaunator.cpp
@@ -2,6 +2,7 @@
#include "delaunator.h"
#include <cstdio>
#include <limits>
+#include <tuple>
using namespace std;
@@ -37,7 +38,45 @@
const double x = (ey * bl - dy * cl) * 0.5 / d;
const double y = (dx * cl - ex * bl) * 0.5 / d;
- return (bl && cl && d && (x * x + y * y)) || max_double;
+ if (bl && cl && d) {
+ return x * x + y * y;
+ } else {
+ return max_double;
+ }
+ }
+
+ double area(
+ const double &px,
+ const double &py,
+ const double &qx,
+ const double &qy,
+ const double &rx,
+ const double &ry
+ ) {
+ return (qy - py) * (rx - qx) - (qx - px) * (ry - qy);
+ }
+
+ tuple<double, double> circumcenter(
+ const double &ax,
+ const double &ay,
+ const double &bx,
+ const double &by,
+ const double &cx,
+ const double &cy
+ ) {
+ const double dx = bx - ax;
+ const double dy = by - ay;
+ const double ex = cx - ax;
+ const double ey = cy - ay;
+
+ const double bl = dx * dx + dy * dy;
+ const double cl = ex * ex + ey * ey;
+ const double d = dx * ey - dy * ex;
+
+ const double x = ax + (ey * bl - dy * cl) * 0.5 / d;
+ const double y = ay + (dx * cl - ex * bl) * 0.5 / d;
+
+ return make_tuple(x, y);
}
}
@@ -87,5 +126,51 @@
min_dist = d;
}
}
- // printf("min_dist=%f i1=%lu", min_dist, i1);
+
+ double min_radius = max_double;
+
+ // find the third point which forms the smallest circumcircle with the first two
+ for (unsigned long int i = 0; i < n; i++)
+ {
+ if (i == i0 || i == i1) continue;
+
+ const double r = circumradius(
+ coords[2 * i0], coords[2 * i0 + 1],
+ coords[2 * i1], coords[2 * i1 + 1],
+ coords[2 * i], coords[2 * i + 1]);
+
+ if (r < min_radius)
+ {
+ i2 = i;
+ min_radius = r;
+ }
+ }
+
+ if (min_radius == max_double) {
+ throw No_triangulation();
+ }
+
+ if (
+ area(
+ coords[2 * i0], coords[2 * i0 + 1],
+ coords[2 * i1], coords[2 * i1 + 1],
+ coords[2 * i2], coords[2 * i2 + 1]
+ ) < 0
+ ) {
+ const double tmp = i1;
+ i1 = i2;
+ i2 = tmp;
+ }
+
+ const double i0x = coords[2 * i0];
+ const double i0y = coords[2 * i0 + 1];
+ const double i1x = coords[2 * i1];
+ const double i1y = coords[2 * i1 + 1];
+ const double i2x = coords[2 * i2];
+ const double i2y = coords[2 * i2 + 1];
+
+ tie(center_x, center_y) = circumcenter(i0x, i0y, i1x, i1y, i2x, i2y);
+ // [cx, cy) = center;
+
+ printf("i0x=%f i0y=%f", i0x, i0y);
};
diff --git a/src/delaunator.h b/src/delaunator.h
index 128f495..bb022bc 100644
--- a/src/delaunator.h
+++ b/src/delaunator.h
@@ -1,9 +1,18 @@
#pragma once
#include <vector>
+#include <exception>
using namespace std;
class Delaunator{
public:
- Delaunator(const vector<double> &coords);
+ Delaunator(const vector<double> &coords);
+ struct No_triangulation : public exception {
+ const char * what () const throw () {
+ return "No delaunay triangulation";
+ };
+ };
+ private:
+ double center_x;
+ double center_y;
};