#![allow(
clippy::cast_precision_loss,
clippy::cast_possible_truncation,
clippy::cast_sign_loss,
reason = "angular vertex-count arithmetic; values are small and non-negative"
)]
#![allow(clippy::float_cmp, reason = "exact degenerate-case guards")]
use alloc::vec::Vec;
use geometry_cs::{CartesianFamily, CoordinateSystem};
use geometry_model::{Point2D, Polygon, Ring};
use geometry_trait::{Point, Polygon as PolygonTrait, Ring as RingTrait};
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum JoinStrategy {
Round {
points_per_circle: usize,
},
Miter,
}
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum PointStrategy {
Circle {
points_per_circle: usize,
},
Square,
}
#[must_use]
pub fn buffer_point<Cs>(
center: &Point2D<f64, Cs>,
distance: f64,
strategy: PointStrategy,
) -> Ring<Point2D<f64, Cs>>
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let cx = center.get::<0>();
let cy = center.get::<1>();
match strategy {
PointStrategy::Circle { points_per_circle } => {
circle_ring(cx, cy, distance, points_per_circle.max(3))
}
PointStrategy::Square => {
let d = distance;
Ring::from_vec(alloc::vec![
Point2D::new(cx - d, cy - d),
Point2D::new(cx + d, cy - d),
Point2D::new(cx + d, cy + d),
Point2D::new(cx - d, cy + d),
Point2D::new(cx - d, cy - d),
])
}
}
}
#[must_use]
pub fn buffer_convex_polygon<Cs>(
polygon: &Polygon<Point2D<f64, Cs>>,
distance: f64,
join: JoinStrategy,
) -> Polygon<Point2D<f64, Cs>>
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let mut verts: Vec<Point2D<f64, Cs>> = distinct_vertices(polygon.exterior());
if verts.len() < 3 {
return Polygon::new(Ring::new());
}
if signed_area_ccw_positive(&verts) < 0.0 {
verts.reverse();
}
let n = verts.len();
let mut boundary: Vec<Point2D<f64, Cs>> = Vec::new();
for i in 0..n {
let prev = verts[(i + n - 1) % n];
let curr = verts[i];
let next = verts[(i + 1) % n];
let n_in = outward_normal(&prev, &curr);
let n_out = outward_normal(&curr, &next);
let p_in = offset(&curr, &n_in, distance);
let p_out = offset(&curr, &n_out, distance);
boundary.push(p_in);
match join {
JoinStrategy::Round { points_per_circle } => {
push_corner_arc(
&mut boundary,
&curr,
&p_in,
&p_out,
distance,
points_per_circle.max(3),
);
}
JoinStrategy::Miter => {
let n_in_ok = n_in.0 != 0.0 || n_in.1 != 0.0;
let n_out_ok = n_out.0 != 0.0 || n_out.1 != 0.0;
let sx = n_in.0 + n_out.0;
let sy = n_in.1 + n_out.1;
let len2 = sx * sx + sy * sy;
if n_in_ok && n_out_ok && len2 > 0.0 {
let scale = 2.0 * distance / len2;
boundary.push(Point2D::new(
curr.get::<0>() + sx * scale,
curr.get::<1>() + sy * scale,
));
}
}
}
boundary.push(p_out);
}
if let Some(first) = boundary.first().copied() {
boundary.push(first);
}
Polygon::new(Ring::from_vec(boundary))
}
fn circle_ring<Cs>(cx: f64, cy: f64, r: f64, segments: usize) -> Ring<Point2D<f64, Cs>>
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let mut pts = Vec::with_capacity(segments + 1);
let step = core::f64::consts::TAU / segments as f64;
for k in 0..segments {
let a = step * k as f64;
pts.push(Point2D::new(cx + r * a.cos(), cy + r * a.sin()));
}
pts.push(pts[0]);
Ring::from_vec(pts)
}
fn distinct_vertices<Cs>(ring: &Ring<Point2D<f64, Cs>>) -> Vec<Point2D<f64, Cs>>
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let mut pts: Vec<Point2D<f64, Cs>> = ring.points().copied().collect();
if pts.len() >= 2 {
let first = pts[0];
let last = pts[pts.len() - 1];
if same(&first, &last) {
pts.pop();
}
}
pts
}
fn signed_area_ccw_positive<Cs>(verts: &[Point2D<f64, Cs>]) -> f64
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let n = verts.len();
let mut acc = 0.0;
for i in 0..n {
let a = &verts[i];
let b = &verts[(i + 1) % n];
acc += a.get::<0>() * b.get::<1>() - b.get::<0>() * a.get::<1>();
}
acc * 0.5
}
fn outward_normal<Cs>(a: &Point2D<f64, Cs>, b: &Point2D<f64, Cs>) -> (f64, f64)
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let dx = b.get::<0>() - a.get::<0>();
let dy = b.get::<1>() - a.get::<1>();
let len = (dx * dx + dy * dy).sqrt();
if len == 0.0 {
return (0.0, 0.0);
}
(dy / len, -dx / len)
}
fn offset<Cs>(p: &Point2D<f64, Cs>, dir: &(f64, f64), distance: f64) -> Point2D<f64, Cs>
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
Point2D::new(
p.get::<0>() + dir.0 * distance,
p.get::<1>() + dir.1 * distance,
)
}
fn push_corner_arc<Cs>(
out: &mut Vec<Point2D<f64, Cs>>,
center: &Point2D<f64, Cs>,
from: &Point2D<f64, Cs>,
to: &Point2D<f64, Cs>,
distance: f64,
points_per_circle: usize,
) where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
let cx = center.get::<0>();
let cy = center.get::<1>();
let a0 = (from.get::<1>() - cy).atan2(from.get::<0>() - cx);
let mut a1 = (to.get::<1>() - cy).atan2(to.get::<0>() - cx);
while a1 < a0 {
a1 += core::f64::consts::TAU;
}
let sweep = a1 - a0;
let steps = ((sweep / core::f64::consts::TAU) * points_per_circle as f64).ceil() as usize;
let steps = steps.max(1);
for k in 1..steps {
let a = a0 + sweep * (k as f64 / steps as f64);
out.push(Point2D::new(
cx + distance * a.cos(),
cy + distance * a.sin(),
));
}
}
fn same<Cs>(a: &Point2D<f64, Cs>, b: &Point2D<f64, Cs>) -> bool
where
Cs: CoordinateSystem<Family = CartesianFamily> + Copy,
{
a.get::<0>() == b.get::<0>() && a.get::<1>() == b.get::<1>()
}
#[cfg(test)]
mod tests {
use super::{JoinStrategy, PointStrategy, buffer_convex_polygon, buffer_point};
use geometry_algorithm::ring_area;
use geometry_cs::Cartesian;
use geometry_model::{Point2D, Polygon, polygon};
use geometry_trait::Polygon as _;
type P = Point2D<f64, Cartesian>;
fn close(a: f64, b: f64, tol: f64) {
assert!((a - b).abs() < tol, "expected {b}, got {a}");
}
#[test]
fn point_circle_area_approximates_pi_r_squared() {
let disc = buffer_point(
&P::new(0.0, 0.0),
2.0,
PointStrategy::Circle {
points_per_circle: 720,
},
);
close(ring_area(&disc).abs(), core::f64::consts::PI * 4.0, 1e-2);
}
#[test]
fn point_square_area() {
let sq = buffer_point(&P::new(0.0, 0.0), 3.0, PointStrategy::Square);
close(ring_area(&sq).abs(), 36.0, 1e-9);
}
#[test]
fn convex_square_round_buffer_area() {
let sq: Polygon<P> = polygon![[(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (0.0, 0.0)]];
let grown = buffer_convex_polygon(
&sq,
1.0,
JoinStrategy::Round {
points_per_circle: 720,
},
);
let expected = 4.0 + 8.0 + core::f64::consts::PI;
close(ring_area(grown.exterior()).abs(), expected, 1e-2);
}
#[test]
fn convex_triangle_round_buffer_grows() {
let tri: Polygon<P> = polygon![[(0.0, 0.0), (4.0, 0.0), (0.0, 3.0), (0.0, 0.0)]];
let base = ring_area(tri.exterior()).abs(); let grown = buffer_convex_polygon(
&tri,
0.5,
JoinStrategy::Round {
points_per_circle: 360,
},
);
assert!(ring_area(grown.exterior()).abs() > base);
}
#[test]
fn buffer_is_winding_independent() {
let ccw: Polygon<P> =
polygon![[(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (0.0, 0.0)]];
let cw: Polygon<P> = polygon![[(0.0, 0.0), (0.0, 2.0), (2.0, 2.0), (2.0, 0.0), (0.0, 0.0)]];
let j = JoinStrategy::Round {
points_per_circle: 720,
};
let expected = 4.0 + 8.0 + core::f64::consts::PI;
let grown_from_counterclockwise =
ring_area(buffer_convex_polygon(&ccw, 1.0, j).exterior()).abs();
let grown_from_clockwise = ring_area(buffer_convex_polygon(&cw, 1.0, j).exterior()).abs();
close(grown_from_counterclockwise, expected, 5e-2);
close(grown_from_clockwise, expected, 5e-2);
}
#[test]
fn miter_square_area_is_16() {
let sq: Polygon<P> = polygon![[(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (0.0, 0.0)]];
let grown = buffer_convex_polygon(&sq, 1.0, JoinStrategy::Miter);
close(ring_area(grown.exterior()).abs(), 16.0, 1e-9);
}
#[test]
fn miter_contains_near_corner_probe() {
use geometry_algorithm::within;
let sq: Polygon<P> = polygon![[(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (0.0, 0.0)]];
let grown = buffer_convex_polygon(&sq, 1.0, JoinStrategy::Miter);
let ang = 22.5_f64.to_radians();
let probe = P::new(2.0 + 0.99 * ang.cos(), 2.0 + 0.99 * ang.sin());
assert!(
within(&probe, &grown),
"buffer must contain points within d"
);
}
#[test]
fn miter_is_superset_of_round_by_area() {
let j_round = JoinStrategy::Round {
points_per_circle: 720,
};
let square: Polygon<P> =
polygon![[(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (0.0, 0.0)]];
let triangle: Polygon<P> = polygon![[(0.0, 0.0), (4.0, 0.0), (0.0, 3.0), (0.0, 0.0)]];
for pg in [square, triangle] {
let m =
ring_area(buffer_convex_polygon(&pg, 1.0, JoinStrategy::Miter).exterior()).abs();
let r = ring_area(buffer_convex_polygon(&pg, 1.0, j_round).exterior()).abs();
assert!(m >= r - 1e-9, "miter {m} must not be below round {r}");
}
}
#[test]
fn miter_is_winding_independent() {
let ccw: Polygon<P> =
polygon![[(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (0.0, 0.0)]];
let cw: Polygon<P> = polygon![[(0.0, 0.0), (0.0, 2.0), (2.0, 2.0), (2.0, 0.0), (0.0, 0.0)]];
close(
ring_area(buffer_convex_polygon(&ccw, 1.0, JoinStrategy::Miter).exterior()).abs(),
16.0,
1e-9,
);
close(
ring_area(buffer_convex_polygon(&cw, 1.0, JoinStrategy::Miter).exterior()).abs(),
16.0,
1e-9,
);
}
}