emh/voronoi
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

voronoi

Browse Source
This commit is contained in:
EvilMuffinHa 2022-06-17 00:27:18 -04:00
parent ac732863d3
commit 07468785b9
3 changed files with 76 additions and 14 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
Cargo.lock generated
View File

@ -139,13 +139,6 @@ dependencies = [
"libc",
]
[[package]]
name = "hilbert"
version = "0.1.0"
dependencies = [
"image",
]
[[package]]
name = "image"
version = "0.23.14"
@ -324,6 +317,13 @@ dependencies = [
"weezl",
]
[[package]]
name = "voronoi"
version = "0.1.0"
dependencies = [
"image",
]
[[package]]
name = "weezl"
version = "0.1.6"

2
Cargo.toml
View File

@ -1,5 +1,5 @@
[package]
name = "hilbert"
name = "voronoi"
version = "0.1.0"
edition = "2021"

74
src/main.rs
View File

@ -27,7 +27,12 @@ impl<const N: usize> Distance<N> for Manhattan<N> {
}
impl<const N: usize> Distance<N> for HilbertPoly {
fn dist(&self, a: [f64; N], b: [f64; N]) -> Result<f64, &'static str> {
if a == b {
return Ok(0.0);
}
// ("#0: {:?} {:?}", a, b);
// 2d Hilbert Convex Polygon, 3d is slightly different as it will be to a side (edge of
// convex set )
@ -57,6 +62,7 @@ impl<const N: usize> Distance<N> for HilbertPoly {
}
// Check if points are both inside the set
let (mut i, mut j, mut c) = (0, self.vertices.len() - 1, false);
loop {
if i >= self.vertices.len() {
@ -94,6 +100,37 @@ impl<const N: usize> Distance<N> for HilbertPoly {
return Err("Second point not in set.");
}
// Check if points are on the boundary
//
let c = edges.clone().into_iter().map(|e| [[e[0][0] - a[0], e[0][1] - a[1]], [e[1][0] - a[0], e[1][1] - a[1]]]).map(|s| {
let ab = f64::sqrt(s[0][0] * s[0][0] + s[0][1] * s[0][1]);
let bc = f64::sqrt(s[1][0] * s[1][0] + s[1][1] * s[1][1]);
let ac = f64::sqrt((s[1][0] - s[0][0]) * (s[1][0] - s[0][0]) + (s[1][1] - s[0][1]) * (s[1][1] - s[0][1]));
// ("#1: boundary check: {} {}", f64::round(ab + bc), f64::round(ac));
if f64::round(ab + bc) != f64::round(ac) {
return true;
}
false
}).all(|x| x);
if !c {
return Err("First point not in set.");
}
let c = edges.clone().into_iter().map(|e| [[e[0][0] - b[0], e[0][1] - b[1]], [e[1][0] - b[0], e[1][1] - b[1]]]).map(|s| {
let ab = f64::sqrt(s[0][0] * s[0][0] + s[0][1] * s[0][1]);
let bc = f64::sqrt(s[1][0] * s[1][0] + s[1][1] * s[1][1]);
let ac = f64::sqrt((s[1][0] - s[0][0]) * (s[1][0] - s[0][0]) + (s[1][1] - s[0][1]) * (s[1][1] - s[0][1]));
if f64::round(ab + bc) != f64::round(ac) {
return true;
}
false
}).all(|x| x);
if !c {
return Err("Second point not in set.");
}
// Intersect the line with each polygon side
let slope = (b[1] - a[1]) / (b[0] - a[0]);
let p_intersect = edges.into_iter().map(|p| {
@ -101,14 +138,17 @@ impl<const N: usize> Distance<N> for HilbertPoly {
let d = edge_slope - slope;
if d == 0.0 {
None
} else if d == f64::INFINITY || d == f64::NEG_INFINITY {
Some([a[0], edge_slope * a[0] - edge_slope * p[0][0] + p[0][1]])
} else {
// ("#2: slopes: {:?} {:?} {} {} {:?}", a, b, slope, edge_slope, d);
Some([((edge_slope * p[0][0] - p[0][1]) - (slope * a[0] - a[1])) / d, ((slope * (edge_slope * p[0][0] - p[0][1])) - (edge_slope * (slope * a[0] - a[1]))) / d])
}
}).filter(|x| x.is_some()).map(|x| x.map_or([0.0, 0.0], |x| x)).collect::<Vec<[f64; 2]>>();
// Find the closest distance
let a_closest = p_intersect.clone().into_iter().map(|p| (p, (p[0] - a[0]) * (p[0] - a[0]) + (p[1] - a[1]) * (p[1] - a[1]))).min_by(|a, b| a.1.partial_cmp(&b.1).unwrap()).unwrap();
let b_closest = p_intersect.into_iter().map(|p| (p, (p[0] - b[0]) * (p[0] - b[0]) + (p[1] - b[1]) * (p[1] - b[1]))).min_by(|a, b| a.1.partial_cmp(&b.1).unwrap()).unwrap();
let b_closest = p_intersect.into_iter().map(|p| (p, (p[0] - b[0]) * (p[0] - b[0]) + (p[1] - b[1]) * (p[1] - b[1]))).filter(|p| p.0 != a_closest.0).min_by(|a, b| a.1.partial_cmp(&b.1).unwrap()).unwrap();
let pb = f64::sqrt((a_closest.0[0] - b[0]) * (a_closest.0[0] - b[0]) + (a_closest.0[1] - b[1]) * (a_closest.0[1] - b[1]));
let pa = f64::sqrt(a_closest.1);
@ -119,23 +159,45 @@ impl<const N: usize> Distance<N> for HilbertPoly {
}
}
fn naive_voronoi_2d(metric: &mut dyn Distance<2>, w: u32, h: u32, points: Vec<(u32, u32)>) {
fn naive_voronoi_2d(metric: &mut dyn Distance<2>, w: u32, h: u32, points: Vec<(u32, u32, Rgb<u8>)>) -> RgbImage {
let mut image: RgbImage = ImageBuffer::new(w, h);
// All points must be in the width / height
for p in points {
for p in points.clone() {
assert!(p.0 <= w);
assert!(p.1 <= h);
}
for i in 0..w {
for j in 0..h {
let c = points.clone().into_iter().map(|p| metric.dist([i as f64, j as f64], [p.0 as f64, p.1 as f64]) /* ("#3: points for dist: ({} {}) ({} {}) {:?}", i, j, p.0, p.1, u)*/);
if c.clone().map(|d| d.is_ok()).all(|x| x) {
let closest = c.clone().map(|d| d.unwrap()).enumerate().min_by(|a, b| a.1.partial_cmp(&b.1).unwrap()).unwrap();
if closest.1 == 0.0 {
*image.get_pixel_mut(i, j) = Rgb([0, 0, 0]);
} else {
*image.get_pixel_mut(i, j) = points[closest.0].2;
}
} else {
*image.get_pixel_mut(i, j) = Rgb([255, 255, 255]);
}
}
}
image
}
fn main() {
let triangle = HilbertPoly { vertices: vec![[-8.0, 0.0], [8.0, 0.0], [0.0, 8.0]] };
println!("{}", triangle.dist([-2.0, 2.0], [2.0, 2.0]).unwrap());
let mut triangle = HilbertPoly { vertices: vec![[300.0, 0.0], [0.0, 600.0], [600.0, 600.0]] };
let mut euclidean: Euclidean<2> = Euclidean();
let mut manhattan: Manhattan<2> = Manhattan();
let points = vec![(360, 340, Rgb([255, 0, 0])), (340, 340, Rgb([0, 255, 0])), (400, 400, Rgb([0, 0, 255]))];
let h = naive_voronoi_2d(&mut triangle, 600, 600, points.clone());
h.save("hilbert.png").unwrap();
let e = naive_voronoi_2d(&mut euclidean, 600, 600, points.clone());
e.save("euclidean.png").unwrap();
let m = naive_voronoi_2d(&mut manhattan, 600, 600, points);
m.save("manhattan.png").unwrap();
}