Skip to main content

ogeom_algo/
seam_columns.rs

1//! Seams stored as plain pcurves, joined.
2//!
3//! An edge along a closed surface's seam that the surface's faces walk
4//! once each way needs both columns of the seam, and a seam
5//! representation holds them so each walk finds its own. Stored instead
6//! as plain pcurves on the one surface (two a period apart, or one or two
7//! on the same column), the lookup by surface finds the first for both
8//! walks, and one of them reads its boundary a period away.
9//! [`join_seam_columns`] puts such an edge right.
10
11use ogeom_core::FastMap;
12
13use ogeom_core::{OgeomResult, Tolerances, ogeom_bail};
14use ogeom_geom::{Curve2d as _, Surface as _};
15use ogeom_math::{Point2, Transform2, Vector2};
16use ogeom_topo::{
17    EdgeData, EdgeRepr, Location, Model, NodeData, Orientation, PCurveId, Shape, ShapeType,
18    SurfaceId, TShapeId, explore_unique,
19};
20
21/// Make a seam of every edge that holds its column on a closed surface as
22/// plain pcurves where the surface's faces walk it twice, once each way,
23/// on two columns a period apart: two faces joined along the surface's
24/// seam, or one face going round it. Returns how many edges were made
25/// seams.
26///
27/// The edge's plain pcurves on the surface, one, or two on one column or a
28/// whole period apart, state where the edge stands in the chart. Each walk
29/// takes the copy of that column, a whole period over or none, whose ends
30/// meet the walk's neighbouring edges. Where the two walks take different
31/// copies, the edge becomes a seam carrying both over one range, each on
32/// the side walked its way, and a lookup by surface finds the right one
33/// for each walk. An edge whose copies the neighbours do not tell apart,
34/// whose walks take one copy, or that the surface's faces walk other than
35/// twice, is left as it is.
36///
37/// Nodes below `shape` that other shapes hold as well are copied first
38/// ([`Model::unshare`]), so the repair reaches no other shape.
39///
40/// # Errors
41///
42/// [`OgeomError::Dangling`](ogeom_core::OgeomError::Dangling) if the shape
43/// or anything it names does not resolve.
44pub fn join_seam_columns(model: &mut Model, shape: &Shape, tol: Tolerances) -> OgeomResult<usize> {
45    model.unshare(shape)?;
46    // Each edge's uses by faces on their surface: the face's boundary ring
47    // as it stores it, and where in the ring the edge stands.
48    let mut uses = FastMap::<(TShapeId, SurfaceId), Vec<(Vec<Shape>, usize)>>::default();
49    for face in explore_unique(model, shape, ShapeType::Face)? {
50        let Some(surface) = model
51            .node(&face)
52            .and_then(|n| n.data().as_face())
53            .map(|d| d.surface)
54        else {
55            ogeom_bail!(Dangling, "face is not in this model");
56        };
57        let stored = face.oriented(Orientation::Forward);
58        for wire in model.ordered_children_of(&stored)? {
59            let ring = model.ordered_children_of(&wire)?;
60            for (at, edge) in ring.iter().enumerate() {
61                uses.entry((edge.node(), surface))
62                    .or_default()
63                    .push((ring.clone(), at));
64            }
65        }
66    }
67    let mut joined = 0;
68    for edge in explore_unique(model, shape, ShapeType::Edge)? {
69        let Some(data) = model.node(&edge).and_then(|n| n.data().as_edge()) else {
70            ogeom_bail!(Dangling, "edge is not in this model");
71        };
72        for stated in stated_columns(model, data, tol)? {
73            let Some([(ring_x, at_x), (ring_y, at_y)]) = uses
74                .get(&(edge.node(), stated.surface))
75                .and_then(|faces| <&[_; 2]>::try_from(&faces[..]).ok())
76            else {
77                continue;
78            };
79            let (walk_x, walk_y) = (ring_x[*at_x].orientation(), ring_y[*at_y].orientation());
80            if walk_x == walk_y {
81                continue;
82            }
83            let (Some(sx), Some(sy)) = (
84                copy_met(model, &stated, ring_x, *at_x, tol)?,
85                copy_met(model, &stated, ring_y, *at_y, tol)?,
86            ) else {
87                continue;
88            };
89            if (sx - sy).magnitude() <= tol.parametric() {
90                continue;
91            }
92            let Some(column) = model.geometry().pcurve(stated.pcurve).cloned() else {
93                ogeom_bail!(Dangling, "an edge names geometry not in this model");
94            };
95            let mut copy_at = |shift: Vector2| -> OgeomResult<PCurveId> {
96                if shift.magnitude() <= tol.parametric() {
97                    return Ok(stated.pcurve);
98                }
99                let moved = column.transformed(&Transform2::translation(shift), tol)?;
100                Ok(model.geometry_mut().add_pcurve(moved))
101            };
102            let (forward, reversed) = if walk_x == Orientation::Forward {
103                (copy_at(sx)?, copy_at(sy)?)
104            } else {
105                (copy_at(sy)?, copy_at(sx)?)
106            };
107            let Some(node) = model.node_mut(&edge) else {
108                ogeom_bail!(Dangling, "edge is not in this model");
109            };
110            let NodeData::Edge(data) = node.data_mut() else {
111                ogeom_bail!(Dangling, "edge node holds no edge data");
112            };
113            // The columns stand where the edge's pcurves stood, or a whole
114            // period over, so its same-parameter claim stands as it was.
115            data.representations[stated.at[0]] = EdgeRepr::Seam {
116                forward,
117                reversed,
118                surface: stated.surface,
119                location: stated.location,
120                range: stated.range,
121            };
122            if let Some(&second) = stated.at.get(1) {
123                data.representations.remove(second);
124            }
125            joined += 1;
126            break;
127        }
128    }
129    Ok(joined)
130}
131
132/// An edge's column on a closed surface as its plain pcurves state it: the
133/// pcurves' places among its representations, the first's pcurve, range
134/// and ends in the chart, and the shifts by whole periods a walk may read
135/// the column at.
136struct StatedColumn {
137    at: Vec<usize>,
138    pcurve: PCurveId,
139    surface: SurfaceId,
140    location: Location,
141    range: (f64, f64),
142    ends: (Point2, Point2),
143    shifts: Vec<Vector2>,
144}
145
146/// The edge's columns on closed surfaces: per surface and placement, one
147/// plain pcurve, or two that stand on one column or a whole period apart
148/// all along, each read over its own range.
149fn stated_columns(
150    model: &Model,
151    data: &EdgeData,
152    tol: Tolerances,
153) -> OgeomResult<Vec<StatedColumn>> {
154    const SAMPLES: u32 = 8;
155    let plain: Vec<_> = data
156        .representations
157        .iter()
158        .enumerate()
159        .filter_map(|(k, r)| match r {
160            EdgeRepr::PCurve {
161                curve,
162                surface,
163                location,
164                range,
165            } => Some((k, *curve, *surface, location, *range)),
166            _ => None,
167        })
168        .collect();
169    let mut stated = Vec::new();
170    for a in &plain {
171        let group: Vec<_> = plain.iter().filter(|p| p.2 == a.2 && p.3 == a.3).collect();
172        // Each group once, from its first member.
173        if group[0].0 != a.0 || group.len() > 2 {
174            continue;
175        }
176        let Some(surface) = model.geometry().surface(a.2) else {
177            ogeom_bail!(Dangling, "an edge names geometry not in this model");
178        };
179        let ((u0, u1), (v0, v1)) = surface.domain();
180        let periods = [
181            (surface.is_periodic_u() || surface.is_closed_u(tol))
182                .then_some(Vector2::new(u1 - u0, 0.0)),
183            (surface.is_periodic_v() || surface.is_closed_v(tol))
184                .then_some(Vector2::new(0.0, v1 - v0)),
185        ];
186        let periods: Vec<Vector2> = periods
187            .into_iter()
188            .flatten()
189            .filter(|p| p.magnitude() > 0.0 && p.magnitude().is_finite())
190            .collect();
191        if periods.is_empty() {
192            continue;
193        }
194        let Some(pa) = model.geometry().pcurve(a.1) else {
195            ogeom_bail!(Dangling, "an edge names geometry not in this model");
196        };
197        if let Some(b) = group.get(1) {
198            let Some(pb) = model.geometry().pcurve(b.1) else {
199                ogeom_bail!(Dangling, "an edge names geometry not in this model");
200            };
201            if !a_whole_period_apart(pa, a.4, pb, b.4, &periods, SAMPLES, tol)? {
202                continue;
203            }
204        }
205        let mut shifts = vec![Vector2::new(0.0, 0.0)];
206        for p in &periods {
207            let more: Vec<Vector2> = shifts.iter().flat_map(|s| [*s + *p, *s - *p]).collect();
208            shifts.extend(more);
209        }
210        stated.push(StatedColumn {
211            at: group.iter().map(|p| p.0).collect(),
212            pcurve: a.1,
213            surface: a.2,
214            location: a.3.clone(),
215            range: a.4,
216            ends: (pa.point_at(a.4.0, tol)?, pa.point_at(a.4.1, tol)?),
217            shifts,
218        });
219    }
220    Ok(stated)
221}
222
223/// Whether `b` stands on `a`'s column or a whole number of `periods` from
224/// it, all along, each read over its own range.
225fn a_whole_period_apart(
226    a: &ogeom_geom::PlanarCurve,
227    ra: (f64, f64),
228    b: &ogeom_geom::PlanarCurve,
229    rb: (f64, f64),
230    periods: &[Vector2],
231    samples: u32,
232    tol: Tolerances,
233) -> OgeomResult<bool> {
234    let mut offset: Option<Vector2> = None;
235    for s in 0..=samples {
236        let f = f64::from(s) / f64::from(samples);
237        let d = b.point_at(rb.0 + (rb.1 - rb.0) * f, tol)?
238            - a.point_at(ra.0 + (ra.1 - ra.0) * f, tol)?;
239        match offset {
240            Some(o) if (d - o).magnitude() > tol.parametric() => return Ok(false),
241            Some(_) => {}
242            None => offset = Some(d),
243        }
244    }
245    let Some(mut rest) = offset else {
246        return Ok(false);
247    };
248    // What is left once every whole period along each closed direction is
249    // taken off.
250    for p in periods {
251        let along = rest.dot(*p) / p.dot(*p);
252        rest -= *p * along.round();
253    }
254    Ok(rest.magnitude() <= tol.parametric())
255}
256
257/// The shift of the edge's column that one walk of it reads: the copy
258/// whose ends meet the edges before and after it in the ring, walked as the
259/// ring walks it; `None` where another copy meets them about as well.
260fn copy_met(
261    model: &Model,
262    stated: &StatedColumn,
263    ring: &[Shape],
264    at: usize,
265    tol: Tolerances,
266) -> OgeomResult<Option<Vector2>> {
267    let n = ring.len();
268    let before = ends_on(model, stated.surface, &ring[(at + n - 1) % n], tol)?;
269    let after = ends_on(model, stated.surface, &ring[(at + 1) % n], tol)?;
270    let nearest = |to: Point2, ends: &[(Point2, Point2)], start: bool| {
271        ends.iter()
272            .map(|e| if start { e.0 } else { e.1 }.distance(to))
273            .fold(f64::INFINITY, f64::min)
274    };
275    let mut misses: Vec<(f64, Vector2)> = stated
276        .shifts
277        .iter()
278        .map(|&shift| {
279            let (mut start, mut end) = (stated.ends.0 + shift, stated.ends.1 + shift);
280            if ring[at].orientation() == Orientation::Reversed {
281                core::mem::swap(&mut start, &mut end);
282            }
283            (
284                nearest(start, &before, false) + nearest(end, &after, true),
285                shift,
286            )
287        })
288        .collect();
289    misses.sort_by(|a, b| a.0.total_cmp(&b.0));
290    // The copy a walk reads meets its neighbours within a fitted trim's
291    // slack; every other stands a period off.
292    match misses[..] {
293        [(near, shift), (far, _), ..]
294            if near.is_finite() && far > 10.0 * near.max(tol.parametric()) =>
295        {
296            Ok(Some(shift))
297        }
298        _ => Ok(None),
299    }
300}
301
302/// Every reading of an edge's ends in a surface's chart, as its occurrence
303/// walks it: each pcurve the edge holds on the surface, each side of a
304/// seam.
305fn ends_on(
306    model: &Model,
307    surface: SurfaceId,
308    edge: &Shape,
309    tol: Tolerances,
310) -> OgeomResult<Vec<(Point2, Point2)>> {
311    let Some(data) = model.node(edge).and_then(|n| n.data().as_edge()) else {
312        ogeom_bail!(Dangling, "edge is not in this model");
313    };
314    let mut ends = Vec::new();
315    for repr in &data.representations {
316        let (ids, range) = match repr {
317            EdgeRepr::PCurve {
318                curve,
319                surface: s,
320                range,
321                ..
322            } if *s == surface => (vec![*curve], *range),
323            EdgeRepr::Seam {
324                forward,
325                reversed,
326                surface: s,
327                range,
328                ..
329            } if *s == surface => (vec![*forward, *reversed], *range),
330            _ => continue,
331        };
332        for id in ids {
333            let Some(pcurve) = model.geometry().pcurve(id) else {
334                ogeom_bail!(Dangling, "an edge names geometry not in this model");
335            };
336            let (a, b) = (
337                pcurve.point_at(range.0, tol)?,
338                pcurve.point_at(range.1, tol)?,
339            );
340            ends.push(if edge.orientation() == Orientation::Reversed {
341                (b, a)
342            } else {
343                (a, b)
344            });
345        }
346    }
347    Ok(ends)
348}