add delaunay related funcs

Cargo.toml

@@ -1,6 +1,6 @@
 [package]
 name = "del-geo"
-version = "0.1.19"
+version = "0.1.20"
 edition = "2021"
 authors = ["Nobuyuki Umetani <[email protected]>"]
 description = "2D/3D geometry utility codes"

src/aabb2.rs

@@ -38,4 +38,19 @@ pub fn signed_distance_aabb<Real>(
     let x_dist = (pos_in.x - x_center).abs() - (max0.x - min0.x) * half;
     let y_dist = (pos_in.y - y_center).abs() - (max0.y - min0.y) * half;
     x_dist.max(y_dist)
+}
+
+pub fn from_vtx2vec<T>(
+    vtx2vec: &[nalgebra::Vector2<T>])
+    -> [T;4]
+    where T: nalgebra::RealField + Copy
+{
+    let mut aabb = [vtx2vec[0][0], vtx2vec[0][1], vtx2vec[0][0], vtx2vec[0][1]];
+    for xy in vtx2vec.iter().skip(1) {
+        if xy[0] < aabb[0] { aabb[0] = xy[0]; }
+        if xy[1] < aabb[1] { aabb[1] = xy[1]; }
+        if xy[0] > aabb[2] { aabb[2] = xy[0]; }
+        if xy[1] > aabb[3] { aabb[3] = xy[1]; }
+    }
+    aabb
 }

src/edge2.rs

@@ -3,10 +3,10 @@
 use num_traits::AsPrimitive;
 
 pub fn culling_intersection_<T>(
-    po_s0: &[T;2],
-    po_e0: &[T;2],
-    po_s1: &[T;2],
-    po_e1: &[T;2]) -> bool
+    po_s0: &[T; 2],
+    po_e0: &[T; 2],
+    po_s1: &[T; 2],
+    po_e1: &[T; 2]) -> bool
     where T: num_traits::Float + 'static + Copy,
           f64: num_traits::AsPrimitive<T>
 {
@@ -28,10 +28,10 @@ pub fn culling_intersection_<T>(
 
 
 pub fn intersection_edge2_<T>(
-    po_s0: &[T;2],
-    po_e0: &[T;2],
-    po_s1: &[T;2],
-    po_e1: &[T;2]) -> Option<(T, T)>
+    po_s0: &[T; 2],
+    po_e0: &[T; 2],
+    po_s1: &[T; 2],
+    po_e1: &[T; 2]) -> Option<(T, T)>
     where T: num_traits::Float + 'static + Copy,
           f64: num_traits::AsPrimitive<T>
 {
@@ -47,14 +47,14 @@ pub fn intersection_edge2_<T>(
 }
 
 pub fn winding_number_<T>(
-    ps: &[T;2],
-    pe: &[T;2],
-    po: &[T;2]) -> T
+    ps: &[T; 2],
+    pe: &[T; 2],
+    po: &[T; 2]) -> T
     where T: num_traits::Float + Copy + 'static,
           f64: AsPrimitive<T>
 {
-    let p0 = crate::vec2::sub_(ps,po);
-    let p1 = crate::vec2::sub_(pe,po);
+    let p0 = crate::vec2::sub_(ps, po);
+    let p1 = crate::vec2::sub_(pe, po);
     let y: T = p1[1] * p0[0] - p1[0] * p0[1];
     let x: T = p0[0] * p1[0] + p0[1] * p1[1];
     y.atan2(x) * std::f64::consts::FRAC_1_PI.as_() * 0.5.as_()
@@ -63,6 +63,16 @@ pub fn winding_number_<T>(
 // -----------------------------
 // below: use nalgebra
 
+pub fn length_squared<T>(
+    p0: &nalgebra::Vector2<T>,
+    p1: &nalgebra::Vector2<T>) -> T
+    where T: std::ops::Sub<Output=T> + std::ops::Mul<Output=T> + std::ops::Add<Output=T> + Copy
+{
+    let x = p0[0] - p1[0];
+    let y = p0[1] - p1[1];
+    x * x + y * y
+}
+
 pub fn intersect_edge2(
     po_s0: &nalgebra::Vector2<f32>,
     po_e0: &nalgebra::Vector2<f32>,

src/tri2.rs

@@ -48,4 +48,39 @@ pub fn dldx_<T>(
     )
 }
 
+//////////
 
+
+pub fn area<T>(
+    v1: &nalgebra::Vector2<T>,
+    v2: &nalgebra::Vector2<T>,
+    v3: &nalgebra::Vector2<T>) -> T
+    where T: num_traits::Float + Copy,
+{
+    ((v2[0] - v1[0]) * (v3[1] - v1[1]) - (v3[0] - v1[0]) * (v2[1] - v1[1])) / (T::one()+T::one())
+}
+
+pub fn circumcenter<T>(
+    p0: &nalgebra::Vector2<T>,
+    p1: &nalgebra::Vector2<T>,
+    p2: &nalgebra::Vector2<T>) -> nalgebra::Vector2::<T>
+where T: num_traits::Float + Copy + std::fmt::Debug,
+{
+    let area = crate::tri2::area(p0, p1, p2);
+    if area.is_zero() { panic!(); }
+    let four = T::one() + T::one() + T::one() + T::one();
+    let sixteen = four * four;
+    let tmp_val = T::one() / (area * area * sixteen);
+
+    let dtmp0 = crate::edge2::length_squared(p1,p2);
+    let dtmp1 = crate::edge2::length_squared(p0,p2);
+    let dtmp2 = crate::edge2::length_squared(p0,p1);
+
+    let etmp0: T = tmp_val * dtmp0 * (dtmp1 + dtmp2 - dtmp0);
+    let etmp1: T = tmp_val * dtmp1 * (dtmp0 + dtmp2 - dtmp1);
+    let etmp2: T = tmp_val * dtmp2 * (dtmp0 + dtmp1 - dtmp2);
+
+    nalgebra::Vector2::<T>::new(
+        etmp0 * p0[0] + etmp1 * p1[0] + etmp2 * p2[0],
+        etmp0 * p0[1] + etmp1 * p1[1] + etmp2 * p2[1])
+}

src/vec2.rs

@@ -13,4 +13,10 @@ pub fn to_na<T>(vtx2xyz: &[T], i_vtx: usize) -> nalgebra::Vector2::<T>
     where T: Copy + nalgebra::RealField
 {
     nalgebra::Vector2::<T>::from_row_slice(&vtx2xyz[i_vtx *2..(i_vtx +1)*2])
+}
+
+pub fn norm_squared<T>(v: &nalgebra::Vector2::<T>) -> T
+where T: std::ops::Mul<Output=T> + std::ops::Add<Output = T> + Copy
+{
+    v[0]*v[0] + v[1]*v[1]
 }