From 1f16aa6a31ae820df13d82e28119aa274989b536 Mon Sep 17 00:00:00 2001 From: Edwin Jakobs Date: Fri, 21 Oct 2022 10:33:24 +0200 Subject: [PATCH] [orx-triangulation] Improve triangulation, add kotlin/js support --- build.gradle | 3 +- orx-jvm/orx-triangulation/build.gradle.kts | 10 - orx-shapes/build.gradle.kts | 4 +- .../README.md | 33 +- orx-triangulation/build.gradle.kts | 72 +++ .../src/commonMain/kotlin/Delaunator.kt | 597 ++++++++++++++++++ .../src/commonMain}/kotlin/Delaunay.kt | 134 ++-- .../kotlin/DelaunayTriangulation.kt | 71 +++ .../src/commonMain/kotlin/DoubleDouble.kt | 320 ++++++++++ .../src/commonMain/kotlin/Predicates.kt | 19 + .../src/commonMain}/kotlin/Voronoi.kt | 278 ++++---- .../src/commonMain/kotlin/VoronoiDiagram.kt | 56 ++ .../src/commonTest/kotlin/TestDelaunay.kt | 68 ++ .../src/commonTest/kotlin/TestVoronoi.kt | 17 + .../commonTest/kotlin/TestVoronoiDiagram.kt | 53 ++ .../src/demo/kotlin/DemoDelaunay01.kt | 20 +- .../src/demo/kotlin/DemoDelaunay02.kt | 12 +- .../src/demo/kotlin/DemoVoronoi01.kt | 16 +- .../src/demo/kotlin/DemoVoronoi02.kt | 33 + .../src/demo/kotlin/DemoVoronoi03.kt | 41 ++ settings.gradle.kts | 2 +- 21 files changed, 1608 insertions(+), 251 deletions(-) delete mode 100644 orx-jvm/orx-triangulation/build.gradle.kts rename {orx-jvm/orx-triangulation => orx-triangulation}/README.md (69%) create mode 100644 orx-triangulation/build.gradle.kts create mode 100644 orx-triangulation/src/commonMain/kotlin/Delaunator.kt rename {orx-jvm/orx-triangulation/src/main => orx-triangulation/src/commonMain}/kotlin/Delaunay.kt (62%) create mode 100644 orx-triangulation/src/commonMain/kotlin/DelaunayTriangulation.kt create mode 100644 orx-triangulation/src/commonMain/kotlin/DoubleDouble.kt create mode 100644 orx-triangulation/src/commonMain/kotlin/Predicates.kt rename {orx-jvm/orx-triangulation/src/main => orx-triangulation/src/commonMain}/kotlin/Voronoi.kt (71%) create mode 100644 orx-triangulation/src/commonMain/kotlin/VoronoiDiagram.kt create mode 100644 orx-triangulation/src/commonTest/kotlin/TestDelaunay.kt create mode 100644 orx-triangulation/src/commonTest/kotlin/TestVoronoi.kt create mode 100644 orx-triangulation/src/commonTest/kotlin/TestVoronoiDiagram.kt rename {orx-jvm/orx-triangulation => orx-triangulation}/src/demo/kotlin/DemoDelaunay01.kt (58%) rename {orx-jvm/orx-triangulation => orx-triangulation}/src/demo/kotlin/DemoDelaunay02.kt (70%) rename {orx-jvm/orx-triangulation => orx-triangulation}/src/demo/kotlin/DemoVoronoi01.kt (61%) create mode 100644 orx-triangulation/src/demo/kotlin/DemoVoronoi02.kt create mode 100644 orx-triangulation/src/demo/kotlin/DemoVoronoi03.kt diff --git a/build.gradle b/build.gradle index 018a9c07..d5cf04ce 100644 --- a/build.gradle +++ b/build.gradle @@ -26,7 +26,8 @@ def multiplatformModules = [ "orx-shapes", "orx-quadtree", "orx-hash-grid", - "orx-depth-camera" + "orx-depth-camera", + "orx-triangulation" ] def doNotPublish = ["openrndr-demos"] diff --git a/orx-jvm/orx-triangulation/build.gradle.kts b/orx-jvm/orx-triangulation/build.gradle.kts deleted file mode 100644 index 3890419e..00000000 --- a/orx-jvm/orx-triangulation/build.gradle.kts +++ /dev/null @@ -1,10 +0,0 @@ -plugins { - org.openrndr.extra.convention.`kotlin-jvm` -} - -dependencies { - api(project(":orx-noise")) - implementation(libs.openrndr.shape) - implementation(libs.openrndr.math) - implementation(libs.delaunator) -} \ No newline at end of file diff --git a/orx-shapes/build.gradle.kts b/orx-shapes/build.gradle.kts index 22ab7e87..89ec5801 100644 --- a/orx-shapes/build.gradle.kts +++ b/orx-shapes/build.gradle.kts @@ -38,7 +38,7 @@ kotlin { @Suppress("UNUSED_VARIABLE") val jvmMain by getting { dependencies { - implementation(project(":orx-jvm:orx-triangulation")) + implementation(project(":orx-triangulation")) } } @@ -58,7 +58,7 @@ kotlin { dependencies { implementation(project(":orx-camera")) implementation(project(":orx-color")) - implementation(project(":orx-jvm:orx-triangulation")) + implementation(project(":orx-triangulation")) } } } diff --git a/orx-jvm/orx-triangulation/README.md b/orx-triangulation/README.md similarity index 69% rename from orx-jvm/orx-triangulation/README.md rename to orx-triangulation/README.md index 765b541f..a65142d4 100644 --- a/orx-jvm/orx-triangulation/README.md +++ b/orx-triangulation/README.md @@ -9,17 +9,17 @@ The functionality comes from a Javascript port of the following libraries: ## Usage -### Delaunay +### DelaunayTriangulation -The entry point is the `Delaunay` class. +The entry point is the `DelaunayTriangulation` class. ```kotlin val points: List - val delaunay = Delaunay.from(points) + val delaunay = DelaunayTriangulation(points) // or - val flatPoints: DoubleArray // (x0, y0, x1, x1, x2, y2) - val delaunay = Delaunay(flatPoints) + + val delaunay = points.delaunayTriangulation() ``` This is how you retrieve the triangulation results: @@ -29,21 +29,19 @@ val triangles: List = delaunay.triangles() val halfedges: List = delaunay.halfedges() val hull: ShapeContour = delaunay.hull() -// Updates the triangulation after the points have been modified in-place. -delaunay.update() ``` ### Voronoi -The bounds specifices where the Voronoi diagram will be clipped. +The bounds specify where the Voronoi diagram will be clipped. ```kotlin val bounds: Rectangle -val delaunay = Delaunay.from(points) -val voronoi = delaunay.voronoi(bounds) +val delaunay = points.delaunayTriangulation() +val voronoi = delaunay.voronoiDiagram(bounds) // or -val voronoi = Voronoi(Delaunay.from(points), bounds) +val voronoi = points.voronoiDiagram(bounds) ``` See [To Infinity and Back Again](https://observablehq.com/@mbostock/to-infinity-and-back-again) for an interactive explanation of Voronoi cell clipping. @@ -51,22 +49,19 @@ See [To Infinity and Back Again](https://observablehq.com/@mbostock/to-infinity- This is how you retrieve th results: ```kotlin -val cells: List = voronoi.cellsPolygons() +val cells: List = voronoi.cellPolygons() val cell: ShapeContour = voronoi.cellPolygon(int) // index -val circumcenters: List = voronoi.circumcenters() +val circumcenters: List = voronoi.circumcenters // Returns true if the cell with the specified index i contains the specified vector -val contaisVector = voronoi.contains(int, Vector2) - -// Updates the Voronoi diagram and underlying triangulation -// after the points have been modified in-place -voronoi.update() +val containsVector = voronoi.contains(int, Vector2) ``` -### Author +### Authors Ricardo Matias / [@ricardomatias](https://github.com/ricardomatias) +Edwin Jakobs / [@edwinRNDR](https://github.com/edwinRNDR) ## Demos ### DemoDelaunay01 diff --git a/orx-triangulation/build.gradle.kts b/orx-triangulation/build.gradle.kts new file mode 100644 index 00000000..64844493 --- /dev/null +++ b/orx-triangulation/build.gradle.kts @@ -0,0 +1,72 @@ +import ScreenshotsHelper.collectScreenshots + +plugins { + org.openrndr.extra.convention.`kotlin-multiplatform` +} + +kotlin { + jvm { + @Suppress("UNUSED_VARIABLE") + val demo by compilations.getting { + // TODO: Move demos to /jvmDemo + defaultSourceSet { + kotlin.srcDir("src/demo/kotlin") + } + collectScreenshots { } + } + compilations.all { + kotlinOptions.jvmTarget = libs.versions.jvmTarget.get() + kotlinOptions.apiVersion = libs.versions.kotlinApi.get() + } + testRuns["test"].executionTask.configure { + useJUnitPlatform() + } + } + js(IR) { + browser() + nodejs() + } + + sourceSets { + @Suppress("UNUSED_VARIABLE") + val commonMain by getting { + dependencies { + api(libs.openrndr.math) + api(libs.openrndr.shape) + } + } + + @Suppress("UNUSED_VARIABLE") + val jvmMain by getting { + } + + @Suppress("UNUSED_VARIABLE") + val jvmDemo by getting { + dependencies { + implementation(project(":orx-shapes")) + implementation(project(":orx-triangulation")) + implementation(project(":orx-noise")) + } + } + + @Suppress("UNUSED_VARIABLE") + val jvmTest by getting { + dependencies { + implementation(kotlin("test-common")) + implementation(kotlin("test-annotations-common")) + implementation(kotlin("test-junit5")) + implementation(libs.kotlin.serialization.json) + runtimeOnly(libs.bundles.jupiter) + implementation(libs.spek.dsl) + implementation(libs.kluent) + } + } + + @Suppress("UNUSED_VARIABLE") + val jsTest by getting { + dependencies { + implementation(kotlin("test-js")) + } + } + } +} \ No newline at end of file diff --git a/orx-triangulation/src/commonMain/kotlin/Delaunator.kt b/orx-triangulation/src/commonMain/kotlin/Delaunator.kt new file mode 100644 index 00000000..efa59df7 --- /dev/null +++ b/orx-triangulation/src/commonMain/kotlin/Delaunator.kt @@ -0,0 +1,597 @@ +package org.openrndr.extra.triangulation + +import kotlin.math.* + +private val EPSILON: Double = 2.0.pow(-52) + +/** + * A Kotlin port of Mapbox's Delaunator incredibly fast JavaScript library for Delaunay triangulation of 2D points. + * + * @description Port of Mapbox's Delaunator (JavaScript) library - https://github.com/mapbox/delaunator + * @property coords flat positions' array - [x0, y0, x1, y1..] + * + * @since f0ed80d - commit + * @author Ricardo Matias + */ +@Suppress("unused") +internal class Delaunator(val coords: DoubleArray) { + private val EDGE_STACK = IntArray(512) + + private var count = coords.size shr 1 + + // arrays that will store the triangulation graph + val maxTriangles = (2 * count - 5).coerceAtLeast(0) + private val _triangles = IntArray(maxTriangles * 3) + private val _halfedges = IntArray(maxTriangles * 3) + + lateinit var triangles: IntArray + lateinit var halfedges: IntArray + + // temporary arrays for tracking the edges of the advancing convex hull + private var hashSize = ceil(sqrt(count * 1.0)).toInt() + private var hullPrev = IntArray(count) // edge to prev edge + private var hullNext = IntArray(count) // edge to next edge + private var hullTri = IntArray(count) // edge to adjacent triangle + private var hullHash = IntArray(hashSize) // angular edge hash + private var hullStart: Int = -1 + + // temporary arrays for sorting points + private var ids = IntArray(count) + private var dists = DoubleArray(count) + + private var cx: Double = Double.NaN + private var cy: Double = Double.NaN + + private var trianglesLen: Int = -1 + + lateinit var hull: IntArray + + init { + update() + } + + fun update() { + if (coords.size <= 2) { + halfedges = IntArray(0) + triangles = IntArray(0) + hull = IntArray(0) + return + } + + + // populate an array of point indices calculate input data bbox + var minX = Double.POSITIVE_INFINITY + var minY = Double.POSITIVE_INFINITY + var maxX = Double.NEGATIVE_INFINITY + var maxY = Double.NEGATIVE_INFINITY + + // points -> points + // minX, minY, maxX, maxY + for (i in 0 until count) { + val x = coords[2 * i] + val y = coords[2 * i + 1] + if (x < minX) minX = x + if (y < minY) minY = y + if (x > maxX) maxX = x + if (y > maxY) maxY = y + + ids[i] = i + } + + val cx = (minX + maxX) / 2 + val cy = (minY + maxY) / 2 + + var minDist = Double.POSITIVE_INFINITY + + var i0: Int = -1 + var i1: Int = -1 + var i2: Int = -1 + + // pick a seed point close to the center + for (i in 0 until count) { + val d = dist(cx, cy, coords[2 * i], coords[2 * i + 1]) + + if (d < minDist) { + i0 = i + minDist = d + } + } + + val i0x = coords[2 * i0] + val i0y = coords[2 * i0 + 1] + + minDist = Double.POSITIVE_INFINITY + + // Find the point closest to the seed + for(i in 0 until count) { + if (i == i0) continue + + val d = dist(i0x, i0y, coords[2 * i], coords[2 * i + 1]) + + if (d < minDist && d > 0) { + i1 = i + minDist = d + } + } + + var i1x = coords[2 * i1] + var i1y = coords[2 * i1 + 1] + + var minRadius = Double.POSITIVE_INFINITY + + // Find the third point which forms the smallest circumcircle with the first two + for (i in 0 until count) { + if(i == i0 || i == i1) continue + + val r = circumradius(i0x, i0y, i1x, i1y, coords[2 * i], coords[2 * i + 1]) + + if(r < minRadius) { + i2 = i + minRadius = r + } + } + + if (minRadius == Double.POSITIVE_INFINITY) { + // order collinear points by dx (or dy if all x are identical) + // and return the list as a hull + for (i in 0 until count) { + val a = (coords[2 * i] - coords[0]) + val b = (coords[2 * i + 1] - coords[1]) + dists[i] = if (a == 0.0) b else a + } + + quicksort(ids, dists, 0, count - 1) + + val nhull = IntArray(count) + var j = 0 + var d0 = Double.NEGATIVE_INFINITY + + for (i in 0 until count) { + val id = ids[i] + if (dists[id] > d0) { + nhull[j++] = id + d0 = dists[id] + } + } + + hull = nhull.copyOf(j) + triangles = IntArray(0) + halfedges = IntArray(0) + + return + } + + var i2x = coords[2 * i2] + var i2y = coords[2 * i2 + 1] + + // swap the order of the seed points for counter-clockwise orientation + if (orient2d(i0x, i0y, i1x, i1y, i2x, i2y) < 0.0) { + val i = i1 + val x = i1x + val y = i1y + i1 = i2 + i1x = i2x + i1y = i2y + i2 = i + i2x = x + i2y = y + } + + + val center = circumcenter(i0x, i0y, i1x, i1y, i2x, i2y) + + this.cx = center[0] + this.cy = center[1] + + for (i in 0 until count) { + dists[i] = dist(coords[2 * i], coords[2 * i + 1], center[0], center[1]) + } + + // sort the points by distance from the seed triangle circumcenter + quicksort(ids, dists, 0, count - 1) + + // set up the seed triangle as the starting hull + hullStart = i0 + var hullSize = 3 + + hullNext[i0] = i1 + hullNext[i1] = i2 + hullNext[i2] = i0 + + hullPrev[i2] = i1 + hullPrev[i0] = i2 + hullPrev[i1] = i0 + + hullTri[i0] = 0 + hullTri[i1] = 1 + hullTri[i2] = 2 + + hullHash.fill(-1) + hullHash[hashKey(i0x, i0y)] = i0 + hullHash[hashKey(i1x, i1y)] = i1 + hullHash[hashKey(i2x, i2y)] = i2 + + trianglesLen = 0 + addTriangle(i0, i1, i2, -1, -1, -1) + + var xp = 0.0 + var yp = 0.0 + + for (k in ids.indices) { + val i = ids[k] + val x = coords[2 * i] + val y = coords[2 * i + 1] + + // skip near-duplicate points + if (k > 0 && abs(x - xp) <= EPSILON && abs(y - yp) <= EPSILON) continue + + xp = x + yp = y + + // skip seed triangle points + if (i == i0 || i == i1 || i == i2) continue + + // find a visible edge on the convex hull using edge hash + var start = 0 + val key = hashKey(x, y) + + for (j in 0 until hashSize) { + start = hullHash[(key + j) % hashSize] + + if (start != -1 && start != hullNext[start]) break + } + + start = hullPrev[start] + + var e = start + var q = hullNext[e] + + while (orient2d(x, y, coords[2 * e], coords[2 * e + 1], coords[2 * q], coords[2 * q + 1]) >= 0) { + e = q + + if (e == start) { + e = -1 + break + } + + q = hullNext[e] + } + + if (e == -1) continue // likely a near-duplicate point skip it + + // add the first triangle from the point + var t = addTriangle(e, i, hullNext[e], -1, -1, hullTri[e]) + + // recursively flip triangles from the point until they satisfy the Delaunay condition + hullTri[i] = legalize(t + 2) + hullTri[e] = t // keep track of boundary triangles on the hull + hullSize++ + + // walk forward through the hull, adding more triangles and flipping recursively + var next = hullNext[e] + q = hullNext[next] + + while (orient2d(x, y, coords[2 * next], coords[2 * next + 1], coords[2 * q], coords[2 * q + 1]) < 0) { + t = addTriangle(next, i, q, hullTri[i], -1, hullTri[next]) + hullTri[i] = legalize(t + 2) + hullNext[next] = next // mark as removed + hullSize-- + + next = q + q = hullNext[next] + } + + // walk backward from the other side, adding more triangles and flipping + if (e == start) { + q = hullPrev[e] + + while (orient2d(x, y, coords[2 * q], coords[2 * q + 1], coords[2 * e], coords[2 * e + 1]) < 0) { + t = addTriangle(q, i, e, -1, hullTri[e], hullTri[q]) + legalize(t + 2) + hullTri[q] = t + hullNext[e] = e // mark as removed + hullSize-- + + e = q + q = hullPrev[e] + } + } + + // update the hull indices + hullStart = e + hullPrev[i] = e + + hullNext[e] = i + hullPrev[next] = i + hullNext[i] = next + + // save the two new edges in the hash table + hullHash[hashKey(x, y)] = i + hullHash[hashKey(coords[2 * e], coords[2 * e + 1])] = e + } + + hull = IntArray(hullSize) + + var e = hullStart + + for (i in 0 until hullSize) { + hull[i] = e + e = hullNext[e] + } + + // trim typed triangle mesh arrays + triangles = _triangles.copyOf(trianglesLen) + halfedges = _halfedges.copyOf(trianglesLen) + } + + private fun legalize(a: Int): Int { + var i = 0 + var na = a + var ar: Int + + // recursion eliminated with a fixed-size stack + while (true) { + val b = _halfedges[na] + + /* if the pair of triangles doesn't satisfy the Delaunay condition + * (p1 is inside the circumcircle of [p0, pl, pr]), flip them, + * then do the same check/flip recursively for the new pair of triangles + * + * pl pl + * /||\ / \ + * al/ || \bl al/ \a + * / || \ / \ + * / a||b \ flip /___ar___\ + * p0\ || /p1 => p0\---bl---/p1 + * \ || / \ / + * ar\ || /br b\ /br + * \||/ \ / + * pr pr + */ + val a0 = na - na % 3 + ar = a0 + (na + 2) % 3 + + if (b == -1) { // convex hull edge + if (i == 0) break + na = EDGE_STACK[--i] + continue + } + + val b0 = b - b % 3 + val al = a0 + (na + 1) % 3 + val bl = b0 + (b + 2) % 3 + + val p0 = _triangles[ar] + val pr = _triangles[na] + val pl = _triangles[al] + val p1 = _triangles[bl] + + val illegal = inCircle( + coords[2 * p0], coords[2 * p0 + 1], + coords[2 * pr], coords[2 * pr + 1], + coords[2 * pl], coords[2 * pl + 1], + coords[2 * p1], coords[2 * p1 + 1]) + + if (illegal) { + _triangles[na] = p1 + _triangles[b] = p0 + + val hbl = _halfedges[bl] + + // edge swapped on the other side of the hull (rare) fix the halfedge reference + if (hbl == -1) { + var e = hullStart + do { + if (hullTri[e] == bl) { + hullTri[e] = na + break + } + e = hullPrev[e] + } while (e != hullStart) + } + link(na, hbl) + link(b, _halfedges[ar]) + link(ar, bl) + + val br = b0 + (b + 1) % 3 + + // don't worry about hitting the cap: it can only happen on extremely degenerate input + if (i < EDGE_STACK.size) { + EDGE_STACK[i++] = br + } + } else { + if (i == 0) break + na = EDGE_STACK[--i] + } + } + + return ar + + } + + private fun link(a:Int, b:Int) { + _halfedges[a] = b + if (b != -1) _halfedges[b] = a + } + + // add a new triangle given vertex indices and adjacent half-edge ids + private fun addTriangle(i0: Int, i1: Int, i2: Int, a: Int, b: Int, c: Int): Int { + val t = trianglesLen + + _triangles[t] = i0 + _triangles[t + 1] = i1 + _triangles[t + 2] = i2 + + link(t, a) + link(t + 1, b) + link(t + 2, c) + + trianglesLen += 3 + + return t + } + + private fun hashKey(x: Double, y: Double): Int { + return (floor(pseudoAngle(x - cx, y - cy) * hashSize) % hashSize).toInt() + } +} + +fun circumradius(ax: Double, ay: Double, + bx: Double, by: Double, + cx: Double, cy: Double): Double { + val dx = bx - ax + val dy = by - ay + val ex = cx - ax + val ey = cy - ay + + val bl = dx * dx + dy * dy + val cl = ex * ex + ey * ey + val d = 0.5 / (dx * ey - dy * ex) + + val x = (ey * bl - dy * cl) * d + val y = (dx * cl - ex * bl) * d + + return x * x + y * y +} + +fun circumcenter(ax: Double, ay: Double, + bx: Double, by: Double, + cx: Double, cy: Double): DoubleArray { + val dx = bx - ax + val dy = by - ay + val ex = cx - ax + val ey = cy - ay + + val bl = dx * dx + dy * dy + val cl = ex * ex + ey * ey + val d = 0.5 / (dx * ey - dy * ex) + + val x = ax + (ey * bl - dy * cl) * d + val y = ay + (dx * cl - ex * bl) * d + + return doubleArrayOf(x, y) +} + +fun quicksort(ids: IntArray, dists: DoubleArray, left: Int, right: Int) { + if (right - left <= 20) { + for (i in (left + 1)..right) { + val temp = ids[i] + val tempDist = dists[temp] + var j = i - 1 + while (j >= left && dists[ids[j]] > tempDist) ids[j + 1] = ids[j--] + ids[j + 1] = temp + } + } else { + val median = (left + right) shr 1 + var i = left + 1 + var j = right + + swap(ids, median, i) + + if (dists[ids[left]] > dists[ids[right]]) swap(ids, left, right) + if (dists[ids[i]] > dists[ids[right]]) swap(ids, i, right) + if (dists[ids[left]] > dists[ids[i]]) swap(ids, left, i) + + val temp = ids[i] + val tempDist = dists[temp] + + while (true) { + do i++ while (dists[ids[i]] < tempDist) + do j-- while (dists[ids[j]] > tempDist) + if (j < i) break + swap(ids, i, j) + } + + ids[left + 1] = ids[j] + ids[j] = temp + + if (right - i + 1 >= j - left) { + quicksort(ids, dists, i, right) + quicksort(ids, dists, left, j - 1) + } else { + quicksort(ids, dists, left, j - 1) + quicksort(ids, dists, i, right) + } + } +} + +private fun swap(arr: IntArray, i: Int, j: Int) { + val tmp = arr[i] + arr[i] = arr[j] + arr[j] = tmp +} + +// monotonically increases with real angle, but doesn't need expensive trigonometry +private fun pseudoAngle(dx: Double, dy: Double): Double { + val p = dx / (abs(dx) + abs(dy)) + val a = if (dy > 0.0) 3.0 - p else 1.0 + p + + return a / 4.0 // [0..1] +} +/* +private fun inCircle(ax: Double, ay: Double, + bx: Double, by: Double, + cx: Double, cy: Double, + px: Double, py: Double): Boolean { + val dx = ax - px + val dy = ay - py + val ex = bx - px + val ey = by - py + val fx = cx - px + val fy = cy - py + + val ap = dx * dx + dy * dy + val bp = ex * ex + ey * ey + val cp = fx * fx + fy * fy + + return dx * (ey * cp - bp * fy) - + dy * (ex * cp - bp * fx) + + ap * (ex * fy - ey * fx) < 0 +}*/ + +private fun inCircle( + ax: Double, ay: Double, + bx: Double, by: Double, + cx: Double, cy: Double, + px: Double, py: Double +): Boolean { + + val dx = twoDiff(ax, px) + val dy = twoDiff(ay, py) + val ex = twoDiff(bx, px) + val ey = twoDiff(by, py) + val fx = twoDiff(cx, px) + val fy = twoDiff(cy, py) + + val ap = ddAddDd(ddMultDd(dx, dx), ddMultDd(dy, dy)) + val bp = ddAddDd(ddMultDd(ex, ex), ddMultDd(ey, ey)) + val cp = ddAddDd(ddMultDd(fx, fx), ddMultDd(fy, fy)) + + val dd = ddAddDd( + ddDiffDd( + ddMultDd(dx, ddDiffDd(ddMultDd(ey, cp), ddMultDd(bp, fy))), + ddMultDd(dy, ddDiffDd(ddMultDd(ex, cp), ddMultDd(bp, fx))) + ), + ddMultDd(ap, ddDiffDd(ddMultDd(ex, fy), ddMultDd(ey, fx))) + ) + // add a small bias here, it seems to help + return (dd[1]) <= 1E-8 +} + + +private fun dist(ax: Double, ay: Double, bx: Double, by: Double): Double { + //val dx = ax - bx + //val dy = ay - by + //return dx * dx + dy * dy + + // double-double implementation but I think it is overkill. + + val dx = twoDiff(ax, bx) + val dy = twoDiff(ay, by) + val dx2 = ddMultDd(dx, dx) + val dy2 = ddMultDd(dy, dy) + val d2 = ddAddDd(dx2, dy2) + + return d2[0] + d2[1] + +} \ No newline at end of file diff --git a/orx-jvm/orx-triangulation/src/main/kotlin/Delaunay.kt b/orx-triangulation/src/commonMain/kotlin/Delaunay.kt similarity index 62% rename from orx-jvm/orx-triangulation/src/main/kotlin/Delaunay.kt rename to orx-triangulation/src/commonMain/kotlin/Delaunay.kt index b1ca53d9..10fec43c 100644 --- a/orx-jvm/orx-triangulation/src/main/kotlin/Delaunay.kt +++ b/orx-triangulation/src/commonMain/kotlin/Delaunay.kt @@ -5,8 +5,10 @@ import org.openrndr.shape.Rectangle import org.openrndr.shape.Triangle import org.openrndr.shape.contour import org.openrndr.shape.contours -import com.github.ricardomatias.Delaunator +import kotlin.js.JsName +import kotlin.math.cos import kotlin.math.pow +import kotlin.math.sin /* ISC License @@ -57,14 +59,14 @@ class Delaunay(val points: DoubleArray) { } } - private var delaunator = Delaunator(points) + private var delaunator: Delaunator = Delaunator(points) val inedges = IntArray(points.size / 2) private val hullIndex = IntArray(points.size / 2) - var halfedges = delaunator.halfedges - var hull = delaunator.hull - var triangles = delaunator.triangles + var halfedges: IntArray = delaunator.halfedges + var hull: IntArray = delaunator.hull + var triangles: IntArray = delaunator.triangles init { init() @@ -75,10 +77,69 @@ class Delaunay(val points: DoubleArray) { init() } + fun neighbors(i:Int) = sequence { + val e0 = inedges.getOrNull(i) ?: return@sequence + if (e0 != -1) { + var e = e0 + var p0 = -1 + + loop@do { + p0 = triangles[e] + yield(p0) + e = if (e % 3 == 2) e - 2 else e + 1 + if (e == -1) { + break@loop + } + + if (triangles[e] != i) { + break@loop + //error("bad triangulation") + } + e = halfedges[e] + + if (e == -1) { + val p = hull[(hullIndex[i] + 1) % hull.size] + if (p != p0) { + yield(p) + } + break@loop + } + } while (e != e0) + } + } + + fun collinear(): Boolean { + for (i in 0 until triangles.size step 3) { + val a = 2 * triangles[i] + val b = 2 * triangles[i + 1] + val c = 2 * triangles[i + 2] + val coords = points + val cross = (coords[c] - coords[a]) * (coords[b + 1] - coords[a + 1]) + - (coords[b] - coords[a]) * (coords[c + 1] - coords[a + 1]) + if (cross > 1e-10) return false; + } + return true + } + private fun jitter(x:Double, y:Double, r:Double): DoubleArray { + return doubleArrayOf(x + sin(x+y) * r, y + cos(x-y)*r) + } fun init() { - halfedges = delaunator.halfedges - hull = delaunator.hull - triangles = delaunator.triangles + + if (hull.size > 2 && collinear()) { + println("warning: triangulation is collinear") + val r = 1E-8 + for (i in 0 until points.size step 2) { + val p = jitter(points[i], points[i+1], r) + points[i] = p[0] + points[i+1] = p[1] + } + + delaunator = Delaunator(points) + halfedges = delaunator.halfedges + hull = delaunator.hull + triangles = delaunator.triangles + + } inedges.fill(-1) hullIndex.fill(-1) @@ -101,52 +162,15 @@ class Delaunay(val points: DoubleArray) { triangles = IntArray(3) { -1 } halfedges = IntArray(3) { -1 } triangles[0] = hull[0] - triangles[1] = hull[1] - triangles[2] = hull[1] inedges[hull[0]] = 1 - if (hull.size == 2) inedges[hull[1]] = 0 + if (hull.size == 2) { + inedges[hull[1]] = 0 + triangles[1] = hull[1] + triangles[2] = hull[1] + } } } - fun triangles(): List { - val list = mutableListOf() - - for (i in triangles.indices step 3 ) { - val t0 = triangles[i] * 2 - val t1 = triangles[i + 1] * 2 - val t2 = triangles[i + 2] * 2 - - val p1 = Vector2(points[t0], points[t0 + 1]) - val p2 = Vector2(points[t1], points[t1 + 1]) - val p3 = Vector2(points[t2], points[t2 + 1]) - - // originally they are defined *counterclockwise* - list.add(Triangle(p3, p2, p1)) - } - - return list - } - - // Inner edges of the delaunay triangulation (without hull) - fun halfedges() = contours { - for (i in halfedges.indices) { - val j = halfedges[i] - - if (j < i) continue - val ti = triangles[i] * 2 - val tj = triangles[j] * 2 - - moveTo(points[ti], points[ti + 1]) - lineTo(points[tj], points[tj + 1]) - } - } - - fun hull() = contour { - for (h in hull) { - moveOrLineTo(points[2 * h], points[2 * h + 1]) - } - close() - } fun find(x: Double, y: Double, i: Int = 0): Int { var i1 = i @@ -178,9 +202,12 @@ class Delaunay(val points: DoubleArray) { c = t } - e = nextHalfedge(e) + e = if (e % 3 == 2) e - 2 else e + 1 - if (triangles[e] != i) break // bad triangulation + if (triangles[e] != i) { + //error("bad triangulation") + break + } // bad triangulation e = halfedges[e] @@ -197,4 +224,5 @@ class Delaunay(val points: DoubleArray) { } fun voronoi(bounds: Rectangle): Voronoi = Voronoi(this, bounds) -} \ No newline at end of file +} + diff --git a/orx-triangulation/src/commonMain/kotlin/DelaunayTriangulation.kt b/orx-triangulation/src/commonMain/kotlin/DelaunayTriangulation.kt new file mode 100644 index 00000000..73a0b314 --- /dev/null +++ b/orx-triangulation/src/commonMain/kotlin/DelaunayTriangulation.kt @@ -0,0 +1,71 @@ +package org.openrndr.extra.triangulation + +import org.openrndr.math.Vector2 +import org.openrndr.shape.Rectangle +import org.openrndr.shape.Triangle +import org.openrndr.shape.contour +import org.openrndr.shape.contours + +/** + * Kotlin/OPENRNDR idiomatic interface to `Delaunay` + */ +class DelaunayTriangulation(val points: List) { + internal val delaunay: Delaunay = Delaunay.from(points) + + fun voronoiDiagram(bounds: Rectangle) = VoronoiDiagram(this, bounds) + + fun neighbors(pointIndex: Int) : Sequence { + return delaunay.neighbors(pointIndex) + } + + fun neighborPoints(pointIndex: Int) : List { + return neighbors(pointIndex).map { points[it] }.toList() + } + + fun triangles(): List { + val list = mutableListOf() + + for (i in delaunay.triangles.indices step 3 ) { + val t0 = delaunay.triangles[i] + val t1 = delaunay.triangles[i + 1] + val t2 = delaunay.triangles[i + 2] + + val p1 = points[t0] + val p2 = points[t1] + val p3 = points[t2] + + // originally they are defined *counterclockwise* + list.add(Triangle(p3, p2, p1)) + } + return list + } + + // Inner edges of the delaunay triangulation (without hull) + fun halfedges() = contours { + for (i in delaunay.halfedges.indices) { + val j = delaunay.halfedges[i] + + if (j < i) continue + val ti = delaunay.triangles[i] + val tj = delaunay.triangles[j] + + moveTo(points[ti]) + lineTo(points[tj]) + } + } + + fun hull() = contour { + for (h in delaunay.hull) { + moveOrLineTo(points[2 * h]) + } + close() + } + + fun nearest(query: Vector2) : Int = delaunay.find(query.x, query.y) + + fun nearestPoint(query: Vector2) : Vector2 = points[nearest(query)] +} + +fun List.delaunayTriangulation() : DelaunayTriangulation { + return DelaunayTriangulation(this) +} \ No newline at end of file diff --git a/orx-triangulation/src/commonMain/kotlin/DoubleDouble.kt b/orx-triangulation/src/commonMain/kotlin/DoubleDouble.kt new file mode 100644 index 00000000..d8650e56 --- /dev/null +++ b/orx-triangulation/src/commonMain/kotlin/DoubleDouble.kt @@ -0,0 +1,320 @@ +package org.openrndr.extra.triangulation + +import kotlin.math.pow + +// original code: https://github.com/FlorisSteenkamp/double-double/ + + +/** + * Returns the difference and exact error of subtracting two floating point + * numbers. + * Uses an EFT (error-free transformation), i.e. `a-b === x+y` exactly. + * The returned result is a non-overlapping expansion (smallest value first!). + * + * * **precondition:** `abs(a) >= abs(b)` - A fast test that can be used is + * `(a > b) === (a > -b)` + * + * See https://people.eecs.berkeley.edu/~jrs/papers/robustr.pdf + */ +internal fun fastTwoDiff(a: Double, b: Double): DoubleArray { + val x = a - b; + val y = (a - x) - b; + + return doubleArrayOf(y, x) +} + +/** + * Returns the sum and exact error of adding two floating point numbers. + * Uses an EFT (error-free transformation), i.e. a+b === x+y exactly. + * The returned sum is a non-overlapping expansion (smallest value first!). + * + * Precondition: abs(a) >= abs(b) - A fast test that can be used is + * (a > b) === (a > -b) + * + * See https://people.eecs.berkeley.edu/~jrs/papers/robustr.pdf + */ +internal fun fastTwoSum(a: Double, b: Double): DoubleArray { + val x = a + b; + + return doubleArrayOf(b - (x - a), x) +} + + +/** + * Truncates a floating point value's significand and returns the result. + * Similar to split, but with the ability to specify the number of bits to keep. + * + * **Theorem 17 (Veltkamp-Dekker)**: Let a be a p-bit floating-point number, where + * p >= 3. Choose a splitting point s such that p/2 <= s <= p-1. Then the + * following algorithm will produce a (p-s)-bit value a_hi and a + * nonoverlapping (s-1)-bit value a_lo such that abs(a_hi) >= abs(a_lo) and + * a = a_hi + a_lo. + * + * * see [Shewchuk](https://people.eecs.berkeley.edu/~jrs/papers/robustr.pdf) + * + * @param a a double + * @param bits the number of significand bits to leave intact + */ +internal fun reduceSignificand( + a: Double, + bits: Int +): Double { + + val s = 53 - bits; + val f = 2.0.pow(s) + 1; + + val c = f * a; + val r = c - (c - a); + + return r; +} + + +/** + * === 2^Math.ceil(p/2) + 1 where p is the # of significand bits in a double === 53. + * @internal + */ +private const val f = 134217729; // 2**27 + 1; + + +/** + * Returns the result of splitting a double into 2 26-bit doubles. + * + * Theorem 17 (Veltkamp-Dekker): Let a be a p-bit floating-point number, where + * p >= 3. Choose a splitting point s such that p/2 <= s <= p-1. Then the + * following algorithm will produce a (p-s)-bit value a_hi and a + * nonoverlapping (s-1)-bit value a_lo such that abs(a_hi) >= abs(a_lo) and + * a = a_hi + a_lo. + * + * see e.g. [Shewchuk](https://people.eecs.berkeley.edu/~jrs/papers/robustr.pdf) + * @param a A double floating point number + */ +private fun split(a: Double): DoubleArray { + val c = f * a; + val a_h = c - (c - a); + val a_l = a - a_h; + + return doubleArrayOf(a_h, a_l) +} + +/** + * Returns the exact result of subtracting b from a. + * + * @param a minuend - a double-double precision floating point number + * @param b subtrahend - a double-double precision floating point number + */ +internal fun twoDiff(a: Double, b: Double): DoubleArray { + val x = a - b; + val bvirt = a - x; + val y = (a - (x + bvirt)) + (bvirt - b); + + return doubleArrayOf(y, x) +} + +/** + * Returns the exact result of multiplying two doubles. + * + * * the resulting array is the reverse of the standard twoSum in the literature. + * + * Theorem 18 (Shewchuk): Let a and b be p-bit floating-point numbers, where + * p >= 6. Then the following algorithm will produce a nonoverlapping expansion + * x + y such that ab = x + y, where x is an approximation to ab and y + * represents the roundoff error in the calculation of x. Furthermore, if + * round-to-even tiebreaking is used, x and y are non-adjacent. + * + * See https://people.eecs.berkeley.edu/~jrs/papers/robustr.pdf + * @param a A double + * @param b Another double + */ +internal fun twoProduct(a: Double, b: Double): DoubleArray { + val x = a * b; + + //const [ah, al] = split(a); + val c = f * a; + val ah = c - (c - a); + val al = a - ah; + //const [bh, bl] = split(b); + val d = f * b; + val bh = d - (d - b); + val bl = b - bh; + + val y = (al * bl) - ((x - (ah * bh)) - (al * bh) - (ah * bl)); + + //const err1 = x - (ah * bh); + //const err2 = err1 - (al * bh); + //const err3 = err2 - (ah * bl); + //const y = (al * bl) - err3; + + return doubleArrayOf(y, x) +} + +internal fun twoSquare(a: Double): DoubleArray { + val x = a * a; + + //const [ah, al] = split(a); + val c = f * a; + val ah = c - (c - a); + val al = a - ah; + + val y = (al * al) - ((x - (ah * ah)) - 2 * (ah * al)); + + return doubleArrayOf(y, x) +} + +/** + * Returns the exact result of adding two doubles. + * + * * the resulting array is the reverse of the standard twoSum in the literature. + * + * Theorem 7 (Knuth): Let a and b be p-bit floating-point numbers. Then the + * following algorithm will produce a nonoverlapping expansion x + y such that + * a + b = x + y, where x is an approximation to a + b and y is the roundoff + * error in the calculation of x. + * + * See https://people.eecs.berkeley.edu/~jrs/papers/robustr.pdf + */ +internal fun twoSum(a: Double, b: Double): DoubleArray { + val x = a + b; + val bv = x - a; + + return doubleArrayOf((a - (x - bv)) + (b - bv), x) +} + +/** + * Returns the result of subtracting the second given double-double-precision + * floating point number from the first. + * + * * relative error bound: 3u^2 + 13u^3, i.e. fl(a-b) = (a-b)(1+ϵ), + * where ϵ <= 3u^2 + 13u^3, u = 0.5 * Number.EPSILON + * * the error bound is not sharp - the worst case that could be found by the + * authors were 2.25u^2 + * + * ALGORITHM 6 of https://hal.archives-ouvertes.fr/hal-01351529v3/document + * @param x a double-double precision floating point number + * @param y another double-double precision floating point number + */ +internal fun ddDiffDd(x: DoubleArray, y: DoubleArray): DoubleArray { + val xl = x[0]; + val xh = x[1]; + val yl = y[0]; + val yh = y[1]; + + //const [sl,sh] = twoSum(xh,yh); + val sh = xh - yh; val _1 = sh - xh; val sl = (xh - (sh - _1)) + (-yh - _1); + //const [tl,th] = twoSum(xl,yl); + val th = xl - yl; val _2 = th - xl; val tl = (xl - (th - _2)) + (-yl - _2); + val c = sl + th; + //const [vl,vh] = fastTwoSum(sh,c) + val vh = sh + c; val vl = c - (vh - sh); + val w = tl + vl + //const [zl,zh] = fastTwoSum(vh,w) + val zh = vh + w; val zl = w - (zh - vh); + + return doubleArrayOf(zl, zh) +} + +/** + * Returns the product of two double-double-precision floating point numbers. + * + * * relative error bound: 7u^2, i.e. fl(a+b) = (a+b)(1+ϵ), + * where ϵ <= 7u^2, u = 0.5 * Number.EPSILON + * the error bound is not sharp - the worst case that could be found by the + * authors were 5u^2 + * + * * ALGORITHM 10 of https://hal.archives-ouvertes.fr/hal-01351529v3/document + * @param x a double-double precision floating point number + * @param y another double-double precision floating point number + */ +internal fun ddMultDd(x: DoubleArray, y: DoubleArray): DoubleArray { + + + //const xl = x[0]; + val xh = x[1]; + //const yl = y[0]; + val yh = y[1]; + + //const [cl1,ch] = twoProduct(xh,yh); + val ch = xh*yh; + val c = f * xh; val ah = c - (c - xh); val al = xh - ah; + val d = f * yh; val bh = d - (d - yh); val bl = yh - bh; + val cl1 = (al*bl) - ((ch - (ah*bh)) - (al*bh) - (ah*bl)); + + //return fastTwoSum(ch,cl1 + (xh*yl + xl*yh)); + val b = cl1 + (xh*y[0] + x[0]*yh); + val xx = ch + b; + + return doubleArrayOf(b - (xx - ch), xx) +} + + +/** + * Returns the result of adding two double-double-precision floating point + * numbers. + * + * * relative error bound: 3u^2 + 13u^3, i.e. fl(a+b) = (a+b)(1+ϵ), + * where ϵ <= 3u^2 + 13u^3, u = 0.5 * Number.EPSILON + * * the error bound is not sharp - the worst case that could be found by the + * authors were 2.25u^2 + * + * ALGORITHM 6 of https://hal.archives-ouvertes.fr/hal-01351529v3/document + * @param x a double-double precision floating point number + * @param y another double-double precision floating point number + */ +internal fun ddAddDd(x: DoubleArray, y: DoubleArray): DoubleArray { + val xl = x[0]; + val xh = x[1]; + val yl = y[0]; + val yh = y[1]; + + //const [sl,sh] = twoSum(xh,yh); + val sh = xh + yh; val _1 = sh - xh; val sl = (xh - (sh - _1)) + (yh - _1); + //val [tl,th] = twoSum(xl,yl); + val th = xl + yl; val _2 = th - xl; val tl = (xl - (th - _2)) + (yl - _2); + val c = sl + th; + //val [vl,vh] = fastTwoSum(sh,c) + val vh = sh + c; val vl = c - (vh - sh); + val w = tl + vl + //val [zl,zh] = fastTwoSum(vh,w) + val zh = vh + w; val zl = w - (zh - vh); + + return doubleArrayOf(zl, zh) +} + + +/** + * Returns the product of a double-double-precision floating point number and a + * double. + * + * * slower than ALGORITHM 8 (one call to fastTwoSum more) but about 2x more + * accurate + * * relative error bound: 1.5u^2 + 4u^3, i.e. fl(a+b) = (a+b)(1+ϵ), + * where ϵ <= 1.5u^2 + 4u^3, u = 0.5 * Number.EPSILON + * * the bound is very sharp + * * probably prefer `ddMultDouble2` due to extra speed + * + * * ALGORITHM 7 of https://hal.archives-ouvertes.fr/hal-01351529v3/document + * @param y a double + * @param x a double-double precision floating point number + */ +internal fun ddMultDouble1(y: Double, x: DoubleArray): DoubleArray { + val xl = x[0]; + val xh = x[1]; + + //val [cl1,ch] = twoProduct(xh,y); + val ch = xh*y; + val c = f * xh; val ah = c - (c - xh); val al = xh - ah; + val d = f * y; val bh = d - (d - y); val bl = y - bh; + val cl1 = (al*bl) - ((ch - (ah*bh)) - (al*bh) - (ah*bl)); + + val cl2 = xl*y; + //val [tl1,th] = fastTwoSum(ch,cl2); + val th = ch + cl2; + val tl1 = cl2 - (th - ch); + + val tl2 = tl1 + cl1; + //val [zl,zh] = fastTwoSum(th,tl2); + val zh = th + tl2; + val zl = tl2 - (zh - th); + + return doubleArrayOf(zl,zh); +} \ No newline at end of file diff --git a/orx-triangulation/src/commonMain/kotlin/Predicates.kt b/orx-triangulation/src/commonMain/kotlin/Predicates.kt new file mode 100644 index 00000000..265c9f83 --- /dev/null +++ b/orx-triangulation/src/commonMain/kotlin/Predicates.kt @@ -0,0 +1,19 @@ +package org.openrndr.extra.triangulation + +internal fun orient2d(bx: Double, by: Double, ax: Double, ay: Double, cx: Double, cy: Double): Double { + // (ax,ay) (bx,by) are swapped such that the sign of the determinant is flipped. which is what Delaunator.kt expects. + + /* + | a b | = | ax - cx ay - cy | + | c d | | bx - cx by - cy | + */ + + val a = ax - cx + val b = ay - cy + val c = bx - cx + val d = by - cy + + val determinant = ddDiffDd(twoProduct(a, d), twoProduct(b, c)) + + return determinant[1] +} \ No newline at end of file diff --git a/orx-jvm/orx-triangulation/src/main/kotlin/Voronoi.kt b/orx-triangulation/src/commonMain/kotlin/Voronoi.kt similarity index 71% rename from orx-jvm/orx-triangulation/src/main/kotlin/Voronoi.kt rename to orx-triangulation/src/commonMain/kotlin/Voronoi.kt index 6f05daf0..b665cca8 100644 --- a/orx-jvm/orx-triangulation/src/main/kotlin/Voronoi.kt +++ b/orx-triangulation/src/commonMain/kotlin/Voronoi.kt @@ -2,8 +2,13 @@ package org.openrndr.extra.triangulation import org.openrndr.math.Vector2 import org.openrndr.shape.Rectangle +import org.openrndr.shape.Shape import org.openrndr.shape.ShapeContour +import org.openrndr.shape.bounds +import kotlin.js.JsName import kotlin.math.abs +import kotlin.math.floor +import kotlin.math.sign /* ISC License @@ -38,7 +43,7 @@ THIS SOFTWARE. class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { private val _circumcenters = DoubleArray(delaunay.points.size * 2) lateinit var circumcenters: DoubleArray - private set + private set val vectors = DoubleArray(delaunay.points.size * 2) @@ -53,9 +58,21 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { fun init() { val points = delaunay.points + + if (delaunay.points.isEmpty()) { + return + } + val triangles = delaunay.triangles val hull = delaunay.hull + if (points.size == 2) { + _circumcenters[0] = points[0] + _circumcenters[1] = points[1] + circumcenters = _circumcenters + return + } + circumcenters = _circumcenters.copyOf(delaunay.triangles.size / 3 * 2) // Compute circumcenters @@ -80,26 +97,20 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { val dy = y2 - y1 val ex = x3 - x1 val ey = y3 - y1 - val bl = dx * dx + dy * dy - val cl = ex * ex + ey * ey val ab = (dx * ey - dy * ex) * 2 - when { - ab == 0.0 -> { - // degenerate case (collinear diagram) - x = (x1 + x3) / 2 - 1e8 * ey - y = (y1 + y3) / 2 + 1e8 * ex - } - abs(ab) < 1e-8 -> { - // almost equal points (degenerate triangle) - x = (x1 + x3) / 2 - y = (y1 + y3) / 2 - } - else -> { - val d = 1 / ab - x = x1 + (ey * bl - dy * cl) * d - y = y1 + (dx * cl - ex * bl) * d - } + if (abs(ab) < 1e-9) { + var a = 1e9 + val r = triangles[0] * 2 + a *= sign((points[r] - x1) * ey - (points[r + 1] - y1) * ex) + x = (x1 + x3) / 2 - a * ey + y = (y1 + y3) / 2 + a * ex + } else { + val d = 1 / ab + val bl = dx * dx + dy * dy + val cl = ex * ex + ey * ey + x = x1 + (ey * bl - dy * cl) * d + y = y1 + (dx * cl - ex * bl) * d } circumcenters[j] = x @@ -142,47 +153,12 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { } - fun cellsPolygons(): List { - val points = delaunay.points - val cells = mutableListOf() - for (i in 0 until (points.size / 2)) { - cellPolygon(i)?.let { - cells.add(it) - } - } - - return cells - } - - fun cellPolygon(i: Int): ShapeContour? { - val points = clip(i) - - if (points == null || points.isEmpty()) return null - - val polygon = mutableListOf(Vector2(points[0], points[1])) - var n = points.size - - while (n > 1 && points[0] == points[n - 2] && points[1] == points[n - 1]) n -= 2 - - for (idx in 2 until n step 2) { - if (points[idx] != points[idx - 2] || points[idx + 1] != points[idx - 1]) { - polygon.add(Vector2(points[idx], points[idx + 1])) - } - } - - return ShapeContour.fromPoints(polygon, true) - } - - fun circumcenters() = circumcenters.toList().windowed(2, 2).map { - Vector2(it[0], it[1]) - } - - fun contains(i: Int, v: Vector2): Boolean { - return contains(i, v.x, v.y) - } private fun cell(i: Int): MutableList? { + + + val inedges = delaunay.inedges val halfedges = delaunay.halfedges val triangles = delaunay.triangles @@ -196,7 +172,7 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { var e = e0 do { - val t = Math.floorDiv(e, 3) // triangle of edge + val t = floor(e / 3.0).toInt() points.add(circumcenters[t * 2]) points.add(circumcenters[t * 2 + 1]) @@ -211,10 +187,44 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { return points } - private fun clip(i: Int): List? { + fun neighbors(i: Int) = sequence { + val ci = clip(i) + if (ci != null) { + for (j in delaunay.neighbors(i)) { + val cj = clip(j) + if (cj != null) { + val li = ci.size + val lj = cj.size + loop@ for (ai in 0 until ci.size step 2) { + for (aj in 0 until cj.size step 2) { + if (ci[ai] == cj[aj] + && ci[ai + 1] == cj[aj + 1] + && ci[(ai + 2) % li] == cj[(aj + lj - 2) % lj] + && ci[(ai + 3) % li] == cj[(aj + lj - 1) % lj] + ) { + yield(j) + break@loop + } + } + } + } + } + } + } + + internal fun clip(i: Int): List? { // degenerate case (1 valid point: return the box) - if (i == 0 && delaunay.hull.size == 1) { - return listOf(bounds.xmax, bounds.ymin, bounds.xmax, bounds.ymax, bounds.xmin, bounds.ymax, bounds.xmin, bounds.ymin) + if (i == 0 && delaunay.points.size == 2) { + return listOf( + bounds.xmax, + bounds.ymin, + bounds.xmax, + bounds.ymax, + bounds.xmin, + bounds.ymax, + bounds.xmin, + bounds.ymin + ) } val points = cell(i) ?: return null @@ -226,7 +236,7 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { val b = !clipVectors[v + 1].isFalsy() return if (a || b) { - this.clipInfinite(i, points, clipVectors[v], clipVectors[v +1], clipVectors[v + 2], clipVectors[v + 3]) + this.clipInfinite(i, points, clipVectors[v], clipVectors[v + 1], clipVectors[v + 2], clipVectors[v + 3]) } else { this.clipFinite(i, points) } @@ -240,40 +250,41 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { vxn: Double, vyn: Double ): List? { - var P: MutableList? = points.mutableCopyOf().also { list -> - // SHAKY - this.project(list[0], list[1], vx0, vy0)?.also { - list.addAll(0, listOf(it.x, it.y)) - } + var P: MutableList? = points.mutableCopyOf() - this.project(list[list.size - 2], list[list.size - 1], vxn, vyn)?.also { - list.addAll(0, listOf(it.x, it.y)) - } - } - - P = clipFinite(i, P!!) + P!! + project(P[0], P[1], vx0, vy0)?.let { p -> P!!.add(0, p[1]); P!!.add(0, p[0]) } + project(P[P.size - 2], P[P.size - 1], vxn, vyn)?.let { p -> P!!.add(p[0]); P!!.add(p[1]) } + P = this.clipFinite(i, P!!) + var n = 0 if (P != null) { - var n = P.size - var c0: Int? + n = P!!.size + var c0 = -1 var c1 = edgeCode(P[n - 2], P[n - 1]) var j = 0 - + var n = P.size while (j < n) { c0 = c1 c1 = edgeCode(P[j], P[j + 1]) - - if ((c0 and c1) != 0) { + if (c0 != 0 && c1 != 0) { j = edge(i, c0, c1, P, j) n = P.size } - j += 2 } - } else if (contains(i, (bounds.xmin + bounds.xmax) / 2, (bounds.ymin + bounds.ymax) / 2)) { - P = mutableListOf(bounds.xmin, bounds.ymin, bounds.xmax, bounds.ymin, bounds.xmax, bounds.ymax, bounds.xmin, bounds.ymax) + } else if (this.contains(i, (bounds.xmin + bounds.xmax) / 2.0, (bounds.ymin + bounds.ymax) / 2.0)) { + P = mutableListOf( + bounds.xmin, + bounds.ymin, + bounds.xmax, + bounds.ymin, + bounds.xmax, + bounds.ymax, + bounds.xmin, + bounds.ymax + ) } - return P } @@ -283,12 +294,12 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { val P = mutableListOf() var x0: Double var y0: Double - var x1= points[n - 2] - var y1= points[n - 1] + var x1 = points[n - 2] + var y1 = points[n - 1] var c0: Int var c1: Int = regionCode(x1, y1) var e0: Int? = null - var e1: Int? = null + var e1: Int? = 0 for (j in 0 until n step 2) { x0 = x1 @@ -314,8 +325,8 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { if (c0 == 0) { S = clipSegment(x0, y0, x1, y1, c0, c1) if (S == null) continue -// sx0 = S[0] -// sy0 = S[1] + sx0 = S[0] + sy0 = S[1] sx1 = S[2] sy1 = S[3] } else { @@ -329,7 +340,7 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { e0 = e1 e1 = this.edgeCode(sx0, sy0) - if (e0.isTruthy() && e1.isTruthy()) this.edge(i, e0!!, e1, P, P.size) + if (e0 != 0 && e1 != 0) this.edge(i, e0!!, e1, P, P.size) P.add(sx0) P.add(sy0) @@ -351,11 +362,19 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { if (e0.isTruthy() && e1.isTruthy()) this.edge(i, e0!!, e1!!, P, P.size); } else if (this.contains(i, (bounds.xmin + bounds.xmax) / 2, (bounds.ymin + bounds.ymax) / 2)) { - return mutableListOf(bounds.xmax, bounds.ymin, bounds.xmax, bounds.ymax, bounds.xmin, bounds.ymax, bounds.xmin, bounds.ymin) + return mutableListOf( + bounds.xmax, + bounds.ymin, + bounds.xmax, + bounds.ymax, + bounds.xmin, + bounds.ymax, + bounds.xmin, + bounds.ymin + ) } else { return null } - return P } @@ -367,7 +386,7 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { var nc0: Int = c0 var nc1: Int = c1 - while(true) { + while (true) { if (nc0 == 0 && nc1 == 0) return doubleArrayOf(nx0, ny0, nx1, ny1) // SHAKY STUFF if ((nc0 and nc1) != 0) return null @@ -408,96 +427,102 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { } private fun regionCode(x: Double, y: Double): Int { - val code = when { + val xcode = when { x < bounds.xmin -> 0b0001 x > bounds.xmax -> 0b0010 else -> 0b0000 } - return code or when { + val ycode = when { y < bounds.ymin -> 0b0100 y > bounds.ymax -> 0b1000 else -> 0b0000 } + return xcode or ycode } private fun contains(i: Int, x: Double, y: Double): Boolean { -// if ((x = +x, x !== x) || (y = +y, y !== y)) return false; + if (x.isNaN() || y.isNaN()) return false return this.delaunay.step(i, x, y) == i; } private fun edge(i: Int, e0: Int, e1: Int, p: MutableList, j: Int): Int { var j = j var e = e0 - while(e != e1) { + loop@ while (e != e1) { var x: Double = Double.NaN var y: Double = Double.NaN - when(e) { + when (e) { 0b0101 -> { // top-left e = 0b0100 - continue + continue@loop } 0b0100 -> { // top e = 0b0110 x = bounds.xmax y = bounds.ymin - break } 0b0110 -> { // top-right e = 0b0010 - continue + continue@loop } 0b0010 -> { // right e = 0b1010 x = bounds.xmax y = bounds.ymax - break } 0b1010 -> { // bottom-right e = 0b1000 - continue + continue@loop } 0b1000 -> { // bottom - e = 0b0001 + e = 0b1001 x = bounds.xmin y = bounds.ymax - break } 0b1001 -> { // bottom-left e = 0b0001 - continue + continue@loop } 0b0001 -> { // left e = 0b0101 x = bounds.xmin y = bounds.ymin - break } } - if ((p[j] != x || p[j + 1] != y) && contains(i, x, y)) { + if (((j < p.size && (p[j] != x)) || ((j + 1) < p.size && p[j + 1] != y)) && contains(i, x, y)) { + require(!x.isNaN()) + require(!y.isNaN()) p.add(j, y) p.add(j, x) j += 2 + } else if (j >= p.size && contains(i, x, y)) { + require(!x.isNaN()) + require(!y.isNaN()) + p.add(x) + p.add(y) + j += 2 } } if (p.size > 4) { var idx = 0 - - while (idx < p.size) { + var n = p.size + while (idx < n) { val j = (idx + 2) % p.size val k = (idx + 4) % p.size - if (p[idx] == p[j] && p[j] == p[k] - || p[idx + 1] == p[j + 1] && p[j + 1] == p[k + 1]) { + if ((p[idx] == p[j] && p[j] == p[k]) + || (p[idx + 1] == p[j + 1] && p[j + 1] == p[k + 1]) + ) { // SHAKY p.removeAt(j) p.removeAt(j) idx -= 2 + n -= 2 } - idx += 2 } } @@ -511,27 +536,25 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { var y = Double.NaN // top - if(vy < 0) { + if (vy < 0) { if (y0 <= bounds.ymin) return null c = (bounds.ymin - y0) / vy - if(c < t) { + if (c < t) { t = c y = bounds.ymin - x = x0 + c * vx + x = x0 + t * vx } - } - // bottom - else if (vy > 0) { + } else if (vy > 0) { // bottom if (y0 >= bounds.ymax) return null c = (bounds.ymax - y0) / vy - if( c < t) { + if (c < t) { t = c y = bounds.ymax - x = x0 + c * vx + x = x0 + t * vx } } // right @@ -545,8 +568,7 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { x = bounds.xmax y = y0 + t * vy } - // left - } else if (vx < 0) { + } else if (vx < 0) { // left if (x0 <= bounds.xmin) return null c = (bounds.xmin - x0) / vx @@ -558,23 +580,23 @@ class Voronoi(val delaunay: Delaunay, val bounds: Rectangle) { } } - if(x.isNaN() || y.isNaN()) return null + if (x.isNaN() || y.isNaN()) return null return Vector2(x, y) } private fun edgeCode(x: Double, y: Double): Int { - val code = when (x) { + val xcode = when (x) { bounds.xmin -> 0b0001 bounds.xmax -> 0b0010 else -> 0b0000 } - - return code or when (y) { + val ycode = when (y) { bounds.ymin -> 0b0100 bounds.ymax -> 0b1000 else -> 0b0000 } + return xcode or ycode } } @@ -592,7 +614,7 @@ private val Rectangle.xmin: Double get() = this.corner.x private val Rectangle.xmax: Double - get() = this.corner.x + width + get() = this.corner.x + width private val Rectangle.ymin: Double get() = this.corner.y @@ -600,4 +622,6 @@ private val Rectangle.ymin: Double private val Rectangle.ymax: Double get() = this.corner.y + height -private fun Double?.isFalsy() = this == null || this == -0.0 || this == 0.0 || isNaN() \ No newline at end of file +private fun Double?.isFalsy() = this == null || this == -0.0 || this == 0.0 || isNaN() + + diff --git a/orx-triangulation/src/commonMain/kotlin/VoronoiDiagram.kt b/orx-triangulation/src/commonMain/kotlin/VoronoiDiagram.kt new file mode 100644 index 00000000..bafb164d --- /dev/null +++ b/orx-triangulation/src/commonMain/kotlin/VoronoiDiagram.kt @@ -0,0 +1,56 @@ +package org.openrndr.extra.triangulation + +import org.openrndr.math.Vector2 +import org.openrndr.shape.Rectangle +import org.openrndr.shape.ShapeContour +import org.openrndr.shape.bounds + +class VoronoiDiagram(val delaunayTriangulation: DelaunayTriangulation, val bounds: Rectangle) { + private val voronoi = Voronoi(delaunayTriangulation.delaunay, bounds) + + val vectors by lazy { + voronoi.vectors.toList().windowed(2, 2).map { + Vector2(it[0], it[1]) + } + } + + val circumcenters by lazy { + voronoi.circumcenters.toList().windowed(2, 2).map { + Vector2(it[0], it[1]) + } + } + + fun cellPolygon(i: Int): ShapeContour { + val points = voronoi.clip(i) + + if (points == null || points.isEmpty()) return ShapeContour.EMPTY + + val polygon = mutableListOf(Vector2(points[0], points[1])) + var n = points.size + + while (n > 1 && points[0] == points[n - 2] && points[1] == points[n - 1]) n -= 2 + + for (idx in 2 until n step 2) { + if (points[idx] != points[idx - 2] || points[idx + 1] != points[idx - 1]) { + polygon.add(Vector2(points[idx], points[idx + 1])) + } + } + return ShapeContour.fromPoints(polygon, true) + } + + fun cellPolygons(): List { + val points = delaunayTriangulation.points + return (points.indices).map { + cellPolygon(it) + } + } + + fun neighbors(cellIndex: Int): Sequence { + return voronoi.neighbors(cellIndex) + } +} + +fun List.voronoiDiagram(bounds: Rectangle = this.bounds): VoronoiDiagram { + val d = this.delaunayTriangulation() + return d.voronoiDiagram(bounds) +} \ No newline at end of file diff --git a/orx-triangulation/src/commonTest/kotlin/TestDelaunay.kt b/orx-triangulation/src/commonTest/kotlin/TestDelaunay.kt new file mode 100644 index 00000000..2cc06909 --- /dev/null +++ b/orx-triangulation/src/commonTest/kotlin/TestDelaunay.kt @@ -0,0 +1,68 @@ +import org.openrndr.extra.triangulation.Delaunay +import org.openrndr.math.Vector2 +import org.openrndr.shape.Circle +import kotlin.test.Test +import kotlin.test.assertEquals +import kotlin.test.assertTrue + +class TestDelaunay { + /** + * Test if an empty triangulation can be made + */ + @Test + fun testEmpty() { + val points = listOf() + val d = Delaunay.from(points) + assertEquals(0, d.triangles.size) + assertEquals(0, d.halfedges.size) + assertEquals(0, d.hull.size) + assertEquals(0, (d.neighbors(0).toList().size)) + } + + /** + * Test if a one point triangulation can be made + */ + @Test + fun testOnePoint() { + val points = listOf(Vector2(100.0, 100.0)) + val d = Delaunay.from(points) + assertEquals(0, (d.neighbors(0).toList().size)) + } + + /** + * Test if a two point triangulation can be made + */ + @Test + fun testTwoPoints() { + val points = listOf(Vector2(100.0, 100.0), Vector2(300.0, 100.0)) + val d = Delaunay.from(points) + println(d.triangles.size) + println("${d.triangles[0]}, ${d.triangles[1]}, ${d.triangles[2]}") + + // this will be one degenerate triangle since we only have 2 points + assertEquals(3, d.triangles.size) + assertEquals(2, d.hull.size) + + assertEquals(1, (d.neighbors(0).toList().size)) + assertEquals(1, (d.neighbors(0).toList().first())) + + assertEquals(1, (d.neighbors(1).toList().size)) + assertEquals(0, (d.neighbors(1).toList().first())) + } + + @Test + fun testThreePointsCollinear() { + val points = listOf(Vector2(100.0, 100.0), Vector2(200.0, 100.0), Vector2(300.0, 100.0)) + val d = Delaunay.from(points) + assertEquals(3, d.triangles.size) + } + + @Test + fun testNeighbors() { + val c = Circle(200.0, 200.0, 150.0).contour.equidistantPositions(20).take(20) + val d = Delaunay.from(c) + for (j in c.indices) { + assertTrue(d.neighbors(j).toList().isNotEmpty()) + } + } +} \ No newline at end of file diff --git a/orx-triangulation/src/commonTest/kotlin/TestVoronoi.kt b/orx-triangulation/src/commonTest/kotlin/TestVoronoi.kt new file mode 100644 index 00000000..59588bd4 --- /dev/null +++ b/orx-triangulation/src/commonTest/kotlin/TestVoronoi.kt @@ -0,0 +1,17 @@ +import org.openrndr.extra.triangulation.Delaunay +import org.openrndr.shape.Circle +import org.openrndr.shape.Rectangle +import kotlin.test.Test +import kotlin.test.assertTrue + +class TestVoronoi { + @Test + fun testNeighbors() { + val c = Circle(200.0, 200.0, 150.0).contour.equidistantPositions(20).take(20) + val d = Delaunay.from(c) + val v = d.voronoi(Rectangle(0.0, 0.0, 400.0, 400.0)) + for (j in c.indices) { + assertTrue(v.neighbors(j).toList().isNotEmpty()) + } + } +} \ No newline at end of file diff --git a/orx-triangulation/src/commonTest/kotlin/TestVoronoiDiagram.kt b/orx-triangulation/src/commonTest/kotlin/TestVoronoiDiagram.kt new file mode 100644 index 00000000..85225d9a --- /dev/null +++ b/orx-triangulation/src/commonTest/kotlin/TestVoronoiDiagram.kt @@ -0,0 +1,53 @@ +import org.openrndr.extra.triangulation.Delaunay +import org.openrndr.extra.triangulation.delaunayTriangulation +import org.openrndr.math.Vector2 +import org.openrndr.shape.Circle +import org.openrndr.shape.Rectangle +import kotlin.test.Test +import kotlin.test.assertEquals +import kotlin.test.assertTrue + +class TestVoronoiDiagram { + @Test + fun testNeighbors() { + val c = Circle(200.0, 200.0, 150.0).contour.equidistantPositions(20).take(20) + val d = Delaunay.from(c) + val v = d.voronoi(Rectangle(0.0, 0.0, 400.0, 400.0)) + for (j in c.indices) { + assertTrue(v.neighbors(j).toList().isNotEmpty()) + } + } + + @Test + fun testEmpty() { + val dt = listOf().delaunayTriangulation() + val v = dt.voronoiDiagram(Rectangle(0.0, 0.0, 400.0, 400.0)) + assertEquals(0, dt.triangles().size) + assertEquals(0, v.cellPolygons().size) + } + + @Test + fun testOnePoint() { + val dt = listOf(Vector2(100.0, 100.0)).delaunayTriangulation() + val v = dt.voronoiDiagram(Rectangle(0.0, 0.0, 400.0, 400.0)) + assertEquals(0, dt.triangles().size) + assertEquals(1, v.cellPolygons().size) + } + + + @Test + fun testTwoPoints() { + val dt = listOf(Vector2(100.0, 100.0), Vector2(300.0, 300.0)).delaunayTriangulation() + val v = dt.voronoiDiagram(Rectangle(0.0, 0.0, 400.0, 400.0)) + assertEquals(1, dt.triangles().size) + assertEquals(2, v.cellPolygons().size) + } + + @Test + fun testThreePointsCollinear() { + val dt = listOf(Vector2(100.0, 100.0), Vector2(200.0, 200.0), Vector2(300.0, 300.0)).delaunayTriangulation() + val v = dt.voronoiDiagram(Rectangle(0.0, 0.0, 400.0, 400.0)) + assertEquals(1, dt.triangles().size) + assertEquals(3, v.cellPolygons().size) + } +} \ No newline at end of file diff --git a/orx-jvm/orx-triangulation/src/demo/kotlin/DemoDelaunay01.kt b/orx-triangulation/src/demo/kotlin/DemoDelaunay01.kt similarity index 58% rename from orx-jvm/orx-triangulation/src/demo/kotlin/DemoDelaunay01.kt rename to orx-triangulation/src/demo/kotlin/DemoDelaunay01.kt index b58da029..c054a79f 100644 --- a/orx-jvm/orx-triangulation/src/demo/kotlin/DemoDelaunay01.kt +++ b/orx-triangulation/src/demo/kotlin/DemoDelaunay01.kt @@ -1,8 +1,7 @@ import org.openrndr.application import org.openrndr.color.ColorRGBa -import org.openrndr.extensions.SingleScreenshot -import org.openrndr.extra.noise.poissonDiskSampling -import org.openrndr.extra.triangulation.Delaunay +import org.openrndr.extra.noise.scatter +import org.openrndr.extra.triangulation.delaunayTriangulation import org.openrndr.math.Vector2 import org.openrndr.shape.Circle @@ -14,28 +13,17 @@ fun main() { title = "Delaunator" } program { - if (System.getProperty("takeScreenshot") == "true") { - extend(SingleScreenshot()) { - this.outputFile = System.getProperty("screenshotPath") - } - } - val circle = Circle(Vector2(400.0), 250.0) + val points = circle.shape.scatter(30.0) - val points = poissonDiskSampling(drawer.bounds, 30.0) - .filter { circle.contains(it) } - - val delaunay = Delaunay.from(points + circle.contour.equidistantPositions(40)) + val delaunay = (points + circle.contour.equidistantPositions(40)).delaunayTriangulation() val triangles = delaunay.triangles().map { it.contour } extend { drawer.clear(ColorRGBa.BLACK) - - for ((i, triangle) in triangles.withIndex()) { drawer.fill = ColorRGBa.PINK.shade(1.0 - i / (triangles.size * 1.2)) drawer.stroke = ColorRGBa.PINK.shade( i / (triangles.size * 1.0) + 0.1) - drawer.contour(triangle) } } diff --git a/orx-jvm/orx-triangulation/src/demo/kotlin/DemoDelaunay02.kt b/orx-triangulation/src/demo/kotlin/DemoDelaunay02.kt similarity index 70% rename from orx-jvm/orx-triangulation/src/demo/kotlin/DemoDelaunay02.kt rename to orx-triangulation/src/demo/kotlin/DemoDelaunay02.kt index e994c33a..28902ba8 100644 --- a/orx-jvm/orx-triangulation/src/demo/kotlin/DemoDelaunay02.kt +++ b/orx-triangulation/src/demo/kotlin/DemoDelaunay02.kt @@ -1,8 +1,7 @@ import org.openrndr.application import org.openrndr.color.ColorRGBa -import org.openrndr.extensions.SingleScreenshot import org.openrndr.extra.noise.poissonDiskSampling -import org.openrndr.extra.triangulation.Delaunay +import org.openrndr.extra.triangulation.delaunayTriangulation import org.openrndr.math.Vector2 import org.openrndr.shape.Rectangle @@ -13,17 +12,10 @@ fun main() { height = 800 } program { - if (System.getProperty("takeScreenshot") == "true") { - extend(SingleScreenshot()) { - this.outputFile = System.getProperty("screenshotPath") - } - } - val frame = Rectangle.fromCenter(Vector2(400.0), 600.0, 600.0) - val points = poissonDiskSampling(frame, 50.0).map { it + frame.corner } - val delaunay = Delaunay.from(points) + val delaunay = points.delaunayTriangulation() val halfedges = delaunay.halfedges() val hull = delaunay.hull() diff --git a/orx-jvm/orx-triangulation/src/demo/kotlin/DemoVoronoi01.kt b/orx-triangulation/src/demo/kotlin/DemoVoronoi01.kt similarity index 61% rename from orx-jvm/orx-triangulation/src/demo/kotlin/DemoVoronoi01.kt rename to orx-triangulation/src/demo/kotlin/DemoVoronoi01.kt index 11603a3f..2b9eebdc 100644 --- a/orx-jvm/orx-triangulation/src/demo/kotlin/DemoVoronoi01.kt +++ b/orx-triangulation/src/demo/kotlin/DemoVoronoi01.kt @@ -1,8 +1,7 @@ import org.openrndr.application import org.openrndr.color.ColorRGBa -import org.openrndr.extensions.SingleScreenshot import org.openrndr.extra.noise.poissonDiskSampling -import org.openrndr.extra.triangulation.Delaunay +import org.openrndr.extra.triangulation.delaunayTriangulation import org.openrndr.math.Vector2 import org.openrndr.shape.Circle import org.openrndr.shape.Rectangle @@ -14,26 +13,19 @@ fun main() { height = 800 } program { - if (System.getProperty("takeScreenshot") == "true") { - extend(SingleScreenshot()) { - this.outputFile = System.getProperty("screenshotPath") - } - } - val circle = Circle(Vector2(400.0), 250.0) val frame = Rectangle.fromCenter(Vector2(400.0), 600.0, 600.0) val points = poissonDiskSampling(drawer.bounds, 30.0) .filter { circle.contains(it) } - val delaunay = Delaunay.from(points + circle.contour.equidistantPositions(40)) - val voronoi = delaunay.voronoi(frame) + val delaunay = (points + circle.contour.equidistantPositions(40)).delaunayTriangulation() + val voronoi = delaunay.voronoiDiagram(frame) - val cells = voronoi.cellsPolygons() + val cells = voronoi.cellPolygons() extend { drawer.clear(ColorRGBa.BLACK) - drawer.fill = null drawer.stroke = ColorRGBa.PINK drawer.contours(cells) diff --git a/orx-triangulation/src/demo/kotlin/DemoVoronoi02.kt b/orx-triangulation/src/demo/kotlin/DemoVoronoi02.kt new file mode 100644 index 00000000..866782bb --- /dev/null +++ b/orx-triangulation/src/demo/kotlin/DemoVoronoi02.kt @@ -0,0 +1,33 @@ +import org.openrndr.application +import org.openrndr.color.ColorRGBa +import org.openrndr.extra.shapes.grid +import org.openrndr.extra.triangulation.delaunayTriangulation +import org.openrndr.shape.Circle + +fun main() { + application { + configure { + width = 720 + height = 720 + } + program { + extend { + val r = drawer.bounds.offsetEdges(-50.0) + val grid = r.grid(8, 8).flatten() + val circles = grid.map { Circle(it.center, it.width / 4.0) } + val points = circles.flatMap { it.contour.equidistantPositions(6) } + drawer.circles(points, 5.0) + val d = points.delaunayTriangulation() + drawer.stroke = ColorRGBa.PINK + drawer.contours(d.halfedges()) + + drawer.stroke = ColorRGBa.YELLOW + drawer.fill = null + drawer.contours(d.voronoiDiagram(drawer.bounds.offsetEdges(-50.0)).cellPolygons()) + + drawer.stroke = ColorRGBa.GRAY + drawer.contours(d.triangles().map { it.contour }) + } + } + } +} \ No newline at end of file diff --git a/orx-triangulation/src/demo/kotlin/DemoVoronoi03.kt b/orx-triangulation/src/demo/kotlin/DemoVoronoi03.kt new file mode 100644 index 00000000..6c652a29 --- /dev/null +++ b/orx-triangulation/src/demo/kotlin/DemoVoronoi03.kt @@ -0,0 +1,41 @@ +import org.openrndr.application +import org.openrndr.color.ColorRGBa +import org.openrndr.extra.shapes.grid +import org.openrndr.extra.triangulation.delaunayTriangulation +import org.openrndr.math.Vector2 +import org.openrndr.math.Vector3 +import org.openrndr.math.transforms.buildTransform +import org.openrndr.shape.Circle + +fun main() { + application { + configure { + width = 750 + height = 1000 + } + program { + extend { + val r = drawer.bounds.offsetEdges(-100.0) + val grid = r.grid(3,6).flatten() + val circles = grid.map { Circle(Vector2.ZERO, 158.975).contour.transform( + buildTransform { + translate(it.center) + rotate(Vector3.UNIT_Z, 0.0) + } + ) } + val points = circles.flatMap { it.contour.equidistantPositions(16).take(16) } + drawer.circles(points, 5.0) + val d = points.delaunayTriangulation() + drawer.stroke = ColorRGBa.PINK + drawer.contours(d.halfedges()) + + drawer.stroke = ColorRGBa.YELLOW + drawer.fill = ColorRGBa.GRAY.opacify(0.5) + drawer.contours(d.voronoiDiagram(drawer.bounds).cellPolygons()) + + drawer.stroke = ColorRGBa.GRAY + drawer.contours(d.triangles().map { it.contour }) + } + } + } +} \ No newline at end of file diff --git a/settings.gradle.kts b/settings.gradle.kts index bea6ab8c..570259c6 100644 --- a/settings.gradle.kts +++ b/settings.gradle.kts @@ -72,7 +72,7 @@ include( "orx-jvm:orx-tensorflow-natives-windows", "orx-timer", "orx-time-operators", - "orx-jvm:orx-triangulation", + "orx-triangulation", "orx-jvm:orx-kinect-common", "orx-jvm:orx-kinect-v1", "orx-jvm:orx-kinect-v1-natives-linux-arm64",