#[allow(unused_imports)]
use super::functions::*;
use super::functions::{bowyer_watson, is_delaunay};
use super::functions::{
circumcircle, circumsphere_3d, clip_polygon_to_bounds, delaunay_to_voronoi,
flip_edge_if_needed, lloyd_relaxation, nearest_site, polygon_area_pt2d, polygon_centroid_pt2d,
};
#[derive(Debug, Clone)]
pub struct DelaunayTriangulation {
pub points: Vec<Point2D>,
pub triangles: Vec<[usize; 3]>,
}
impl DelaunayTriangulation {
pub fn from_points(pts: &[[f64; 2]]) -> Self {
let points: Vec<Point2D> = pts.iter().map(|&p| Point2D::from_array(p)).collect();
let raw = bowyer_watson(pts);
let triangles: Vec<[usize; 3]> = raw.iter().map(|t| [t.a, t.b, t.c]).collect();
Self { points, triangles }
}
pub fn circumscribed_circle(&self, tri: [usize; 3]) -> ([f64; 2], f64) {
let pa = self.points[tri[0]].to_array();
let pb = self.points[tri[1]].to_array();
let pc = self.points[tri[2]].to_array();
match circumcircle(pa, pb, pc) {
Some((center, r2)) => (center, r2.sqrt()),
None => ([0.0, 0.0], 0.0),
}
}
pub fn is_delaunay(&self) -> bool {
let pts: Vec<[f64; 2]> = self.points.iter().map(|p| p.to_array()).collect();
let tris: Vec<DelaunayTriangle> = self
.triangles
.iter()
.map(|&t| DelaunayTriangle::new(t[0], t[1], t[2]))
.collect();
is_delaunay(&pts, &tris)
}
pub fn flip_edge(&mut self) {
let pts: Vec<[f64; 2]> = self.points.iter().map(|p| p.to_array()).collect();
let mut tris: Vec<DelaunayTriangle> = self
.triangles
.iter()
.map(|&t| DelaunayTriangle::new(t[0], t[1], t[2]))
.collect();
let n = tris.len();
for i in 0..n {
for j in (i + 1)..n {
if j < tris.len() {
flip_edge_if_needed(&mut tris, &pts, i, j);
}
}
}
self.triangles = tris.iter().map(|t| [t.a, t.b, t.c]).collect();
}
pub fn dual_voronoi(&self) -> VoronoiDiagram {
let pts: Vec<[f64; 2]> = self.points.iter().map(|p| p.to_array()).collect();
if pts.is_empty() {
return VoronoiDiagram {
sites: vec![],
cells: vec![],
bbox: [0.0, 1.0, 0.0, 1.0],
};
}
let (mut xmin, mut xmax, mut ymin, mut ymax) = (f64::MAX, f64::MIN, f64::MAX, f64::MIN);
for &p in &pts {
if p[0] < xmin {
xmin = p[0];
}
if p[0] > xmax {
xmax = p[0];
}
if p[1] < ymin {
ymin = p[1];
}
if p[1] > ymax {
ymax = p[1];
}
}
let margin_x = (xmax - xmin).max(1.0) * 0.2;
let margin_y = (ymax - ymin).max(1.0) * 0.2;
let bbox = [
xmin - margin_x,
xmax + margin_x,
ymin - margin_y,
ymax + margin_y,
];
VoronoiDiagram::from_sites(&pts, bbox)
}
}
#[derive(Debug, Clone)]
pub struct VoronoiDiagram {
pub sites: Vec<VoronoiSite>,
pub cells: Vec<VoronoiCell>,
pub bbox: [f64; 4],
}
impl VoronoiDiagram {
pub fn from_sites(pts: &[[f64; 2]], bbox: [f64; 4]) -> Self {
let sites: Vec<VoronoiSite> = pts
.iter()
.enumerate()
.map(|(i, &p)| VoronoiSite::new(Point2D::from_array(p), i))
.collect();
if pts.is_empty() {
return Self {
sites,
cells: vec![],
bbox,
};
}
let triangles = bowyer_watson(pts);
let legacy_cells = delaunay_to_voronoi(pts, &triangles);
let cells: Vec<VoronoiCell> = legacy_cells
.iter()
.enumerate()
.map(|(i, lc)| {
let clipped = clip_polygon_to_bounds(&lc.vertices, bbox);
let verts: Vec<Point2D> = clipped.iter().map(|&v| Point2D::from_array(v)).collect();
VoronoiCell::new(i, verts)
})
.collect();
Self { sites, cells, bbox }
}
pub fn nearest_site(&self, p: [f64; 2]) -> usize {
if self.sites.is_empty() {
return 0;
}
let pts: Vec<[f64; 2]> = self.sites.iter().map(|s| s.pos.to_array()).collect();
nearest_site(p, &pts)
}
pub fn cell_area(&self, i: usize) -> f64 {
if i < self.cells.len() {
self.cells[i].area
} else {
0.0
}
}
pub fn centroid(&self, i: usize) -> [f64; 2] {
if i < self.cells.len() {
let c = self.cells[i].centroid();
[c.x, c.y]
} else {
[0.0, 0.0]
}
}
pub fn lloyd_relaxation(&self, n_iter: usize) -> Self {
let pts_arr: Vec<[f64; 2]> = self.sites.iter().map(|s| s.pos.to_array()).collect();
let relaxed = lloyd_relaxation(&pts_arr, self.bbox, n_iter);
Self::from_sites(&relaxed, self.bbox)
}
}
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct Point2D {
pub x: f64,
pub y: f64,
}
impl Point2D {
pub fn new(x: f64, y: f64) -> Self {
Self { x, y }
}
#[inline]
pub fn dist2(&self, other: &Point2D) -> f64 {
(self.x - other.x).powi(2) + (self.y - other.y).powi(2)
}
#[inline]
pub fn dist(&self, other: &Point2D) -> f64 {
self.dist2(other).sqrt()
}
#[inline]
pub fn to_array(&self) -> [f64; 2] {
[self.x, self.y]
}
#[inline]
pub fn from_array(a: [f64; 2]) -> Self {
Self { x: a[0], y: a[1] }
}
}
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub struct DelaunayTriangle {
pub a: usize,
pub b: usize,
pub c: usize,
}
impl DelaunayTriangle {
pub fn new(a: usize, b: usize, c: usize) -> Self {
Self { a, b, c }
}
pub fn indices(&self) -> [usize; 3] {
[self.a, self.b, self.c]
}
pub fn contains_super_vertex(&self, offset: usize) -> bool {
self.a >= offset || self.b >= offset || self.c >= offset
}
}
#[derive(Debug, Clone)]
pub struct VoronoiCell {
pub site_index: usize,
pub vertices: Vec<Point2D>,
pub area: f64,
}
impl VoronoiCell {
pub fn new(site_index: usize, vertices: Vec<Point2D>) -> Self {
let area = polygon_area_pt2d(&vertices);
Self {
site_index,
vertices,
area,
}
}
pub fn centroid(&self) -> Point2D {
polygon_centroid_pt2d(&self.vertices)
}
pub fn perimeter(&self) -> f64 {
let n = self.vertices.len();
if n < 2 {
return 0.0;
}
(0..n)
.map(|i| self.vertices[i].dist(&self.vertices[(i + 1) % n]))
.sum()
}
}
#[derive(Debug, Clone)]
pub struct LegacyVoronoiCell {
pub site: [f64; 2],
pub vertices: Vec<[f64; 2]>,
}
impl LegacyVoronoiCell {
pub fn new(site: [f64; 2], vertices: Vec<[f64; 2]>) -> Self {
Self { site, vertices }
}
}
#[allow(dead_code)]
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub struct DelaunayTetrahedron {
pub v: [usize; 4],
}
impl DelaunayTetrahedron {
pub fn new(a: usize, b: usize, c: usize, d: usize) -> Self {
Self { v: [a, b, c, d] }
}
pub fn contains_super(&self, offset: usize) -> bool {
self.v.iter().any(|&i| i >= offset)
}
pub fn circumsphere(&self, pts: &[[f64; 3]]) -> Option<([f64; 3], f64)> {
circumsphere_3d(
pts[self.v[0]],
pts[self.v[1]],
pts[self.v[2]],
pts[self.v[3]],
)
}
}
#[allow(dead_code)]
#[derive(Debug, Clone, Copy)]
pub struct WeightedSite {
pub position: [f64; 2],
pub weight: f64,
}
impl WeightedSite {
pub fn new(position: [f64; 2], weight: f64) -> Self {
Self { position, weight }
}
pub fn power_distance(&self, q: [f64; 2]) -> f64 {
let dx = q[0] - self.position[0];
let dy = q[1] - self.position[1];
dx * dx + dy * dy - self.weight
}
}
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct VoronoiSite {
pub pos: Point2D,
pub index: usize,
}
impl VoronoiSite {
pub fn new(pos: Point2D, index: usize) -> Self {
Self { pos, index }
}
}