geop_lib_booleans/
section.rs

1//! Section curves: where two bodies' boundaries cross, or a body's crosses
2//! a plane, as a wire of edges standing on their own.
3//!
4//! Nothing new finds them. A copy of each body is imprinted onto the other
5//! exactly as a boolean imprints two solids ([`remesh`]), after which every
6//! curve where their boundaries cross is an edge both copies share. Those
7//! edges, with the vertices they run between, become the wire; the copies
8//! go, and the bodies are left as they are.
9
10use std::collections::HashMap;
11
12use geop_core_geometry::{aabb::union_aabb, nurb_surface::NurbSurface3D};
13use geop_core_math::{
14    geop_error::{GeopError, GeopResult, WithContext},
15    primitives::CoordinateSystem,
16    scalars::Scalar,
17    vector::Vector2,
18    with_context,
19};
20use geop_core_topology::{
21    Body, VertexId, WireId,
22    build::{BodySpec, EdgeSpec},
23};
24use geop_ops_model::{BodyNames, NamedModel, Namer};
25
26use crate::{
27    imprint::shared_edges,
28    remesh::{RemeshParams, remesh},
29};
30
31/// The curves where the boundaries of the bodies `a` and `b` cross, as a
32/// wire, both bodies left as they are. Where a face of one lies on a face
33/// of the other, the outline of what they share is among the curves.
34///
35/// Each body is imprinted as a copy, the copy of the entity `X` named
36/// `N(X)`; a curve, and a vertex it ends at, is then named as a boolean of
37/// the copies names it (see `naming`), and the wire's edges and vertices
38/// keep those names.
39///
40/// Fails, leaving the model as it was, if the bodies' boundaries do not
41/// meet: there is nothing to make a wire of.
42pub fn section(
43    model: &mut NamedModel,
44    namer: &Namer,
45    a: Body,
46    b: Body,
47    params: RemeshParams,
48) -> GeopResult<WireId> {
49    let before = model.clone();
50    section_of(model, namer, a, b, params).inspect_err(|_| *model = before)
51}
52
53/// [`section`], leaving what it got to in the model where it fails.
54fn section_of(
55    model: &mut NamedModel,
56    namer: &Namer,
57    a: Body,
58    b: Body,
59    params: RemeshParams,
60) -> GeopResult<WireId> {
61    let ctx = with_context!("section(name={}, a={a}, b={b})", namer.root());
62    let copies = [a, b]
63        .map(|body| {
64            model
65                .copy_body(body, |name| namer.name(&[name]))
66                .and_then(|built| built.body())
67        })
68        .into_iter()
69        .collect::<GeopResult<Vec<Body>>>()
70        .with_context(&ctx)?;
71    let [copy_a, copy_b] = copies[..] else {
72        unreachable!("two bodies copied")
73    };
74    remesh(model, namer, copy_a, copy_b, params).with_context(&ctx)?;
75    let curves = shared_edges(model, copy_a, copy_b).with_context(&ctx)?;
76    if curves.is_empty() {
77        return Err(ctx(GeopError::new(
78            "section: the bodies' boundaries do not meet",
79        )));
80    }
81
82    // The wire: every curve, between the vertices it ends at.
83    let topology = model.topology();
84    let mut spec = BodySpec::new(false);
85    let mut names = BodyNames::default();
86    let mut index: HashMap<VertexId, usize> = HashMap::new();
87    let mut vertex = |spec: &mut BodySpec, names: &mut BodyNames, v: VertexId| {
88        Ok::<_, GeopError>(match index.get(&v) {
89            Some(&i) => i,
90            None => {
91                spec.vertices.push(topology.get_vertex(v)?.point);
92                names.vertices.push(model.name(v)?.to_string());
93                index.insert(v, spec.vertices.len() - 1);
94                spec.vertices.len() - 1
95            }
96        })
97    };
98    for &edge in &curves {
99        let e = topology.get_edge(edge)?;
100        let start = vertex(&mut spec, &mut names, e.start_vertex)?;
101        let end = vertex(&mut spec, &mut names, e.end_vertex)?;
102        spec.edges.push(EdgeSpec {
103            curve: e.curve.clone(),
104            start,
105            end,
106        });
107        names.edges.push(model.name(edge)?.to_string());
108    }
109    model.delete_bodies(&[copy_a, copy_b]).with_context(&ctx)?;
110    let built = model.build_body(spec, names).with_context(&ctx)?;
111    built
112        .wire
113        .ok_or_else(|| ctx(GeopError::new("section: the curves built no wire")))
114}
115
116/// The curves where the boundary of the body `body` crosses the plane
117/// `plane` (through its origin, along its `u` and `v`), as a wire, the body
118/// left as it is: [`section`] of the body and a square of the plane reaching
119/// past the whole body. The square's corners and sides are named
120/// `N(plane,v0)` ... and `N(plane,e0)` ..., its face `N(plane)`. Fails, as
121/// [`section`] does, leaving the model as it was.
122pub fn section_by_plane(
123    model: &mut NamedModel,
124    namer: &Namer,
125    body: Body,
126    plane: &CoordinateSystem,
127    params: RemeshParams,
128) -> GeopResult<WireId> {
129    let ctx = with_context!("section_by_plane(name={}, body={body})", namer.root());
130    let topology = model.topology();
131    let faces = topology.body_faces(body).with_context(&ctx)?;
132    let bounds = union_aabb(faces.iter().map(|&f| topology.faces[&f].surface.aabb()))
133        .ok_or_else(|| ctx(GeopError::new("section: the body has no face")))?;
134    // The square, centred where the body's box's centre falls on the plane,
135    // reaching past the box's diagonal each way: a free choice, as long as
136    // the square's own sides are clear of the body.
137    let center = geop_core_math::vector::Vector3::from_array(bounds.map(|x| x.midpoint()));
138    let diagonal = bounds
139        .iter()
140        .map(|x| x.width().mul(x.width()))
141        .fold(Scalar::ZERO, |sum, w| sum.add(w))
142        .sqrt()?;
143    let reach = Scalar::from_f64((diagonal.upper().to_f64() + 1.0).ceil());
144    let middle = plane.to_uv(&center);
145    let (ru, rv) = (
146        reach.div(plane.u().norm())?.sharpen(),
147        reach.div(plane.v().norm())?.sharpen(),
148    );
149    let corner = |su: Scalar, sv: Scalar| {
150        plane.uv_to_xyz(&Vector2::from_array([
151            middle[0].sharpen().add(su.mul(ru)),
152            middle[1].sharpen().add(sv.mul(rv)),
153        ]))
154    };
155    let (minus, plus) = (Scalar::ONE.neg(), Scalar::ONE);
156    let square = NurbSurface3D::bilinear(
157        corner(minus, minus),
158        corner(plus, minus),
159        corner(plus, plus),
160        corner(minus, plus),
161    )
162    .with_context(&ctx)?;
163    let spec = BodySpec::sheet_on(square).with_context(&ctx)?;
164    let names = BodyNames {
165        vertices: (0..4)
166            .map(|k| namer.name(&["plane", &format!("v{k}")]))
167            .collect(),
168        edges: (0..4)
169            .map(|k| namer.name(&["plane", &format!("e{k}")]))
170            .collect(),
171        faces: vec![namer.name(&["plane"])],
172        solid: None,
173    };
174    let sheet = model
175        .build_body(spec, names)
176        .and_then(|built| built.body())
177        .with_context(&ctx)?;
178    let wire = section(model, namer, body, sheet, params).with_context(&ctx);
179    model.delete_bodies(&[sheet]).with_context(&ctx)?;
180    wire
181}
182
183#[cfg(test)]
184mod tests {
185    use super::{section, section_by_plane};
186    use crate::remesh::RemeshParams;
187    use geop_core_math::{primitives::CoordinateSystem, scalars::Scalar};
188    use geop_core_topology::{Body, Sense};
189    use geop_lib_sweep::shapes::{Axis, block, revolved_cylinder_along_axis, v3};
190    use geop_ops_model::{NamedModel, Namer, testing::assert_valid};
191
192    /// Two unit cubes, the second moved by half along each axis, cross in
193    /// one loop: six edges, three across each cube's faces, joined at the
194    /// six points where an edge of one pierces a face of the other. Both
195    /// cubes stay as they were.
196    #[test]
197    fn two_overlapping_cubes_cross_in_one_loop() {
198        let mut model = NamedModel::new();
199        let a = block(&mut model, "a", [0.0; 3], [1.0; 3]);
200        let b = block(&mut model, "b", [0.5; 3], [1.5; 3]);
201        let (faces, edges) = (model.topology().faces.len(), model.topology().edges.len());
202        let namer = Namer::new("section", "s").unwrap();
203        let wire = section(
204            &mut model,
205            &namer,
206            a.into(),
207            b.into(),
208            RemeshParams::default(),
209        )
210        .unwrap();
211        assert_valid(&model);
212        assert_eq!(
213            model.topology().faces.len(),
214            faces,
215            "the cubes as they were"
216        );
217        assert_eq!(model.topology().edges.len(), edges + 6);
218        let chains = model.topology().wire_chains(wire).unwrap();
219        assert_eq!(chains.len(), 1);
220        assert_eq!(chains[0].len(), 6);
221        let topology = model.topology();
222        let end = |&(e, sense): &(_, Sense)| {
223            let e = topology.get_edge(e).unwrap();
224            match sense {
225                Sense::Forward => (e.start_vertex, e.end_vertex),
226                Sense::Reversed => (e.end_vertex, e.start_vertex),
227            }
228        };
229        let ends: Vec<_> = chains[0].iter().map(end).collect();
230        for k in 0..6 {
231            assert_eq!(ends[k].1, ends[(k + 1) % 6].0, "a closed chain");
232        }
233    }
234
235    /// A plane across a cylinder's axis crosses its side in a circle, one
236    /// edge per quarter of the side.
237    #[test]
238    fn a_plane_across_a_cylinder_cuts_a_circle() {
239        let mut model = NamedModel::new();
240        let cylinder = revolved_cylinder_along_axis(
241            &mut model,
242            "c",
243            v3(0.0, 0.0, 0.0),
244            Scalar::from_f64(0.5),
245            Scalar::from_f64(2.0),
246            Axis::Z,
247        )
248        .unwrap();
249        let namer = Namer::new("section", "s").unwrap();
250        let plane = CoordinateSystem::world_at(v3(0.0, 0.0, 0.7));
251        let wire = section_by_plane(
252            &mut model,
253            &namer,
254            Body::Solid(cylinder),
255            &plane,
256            RemeshParams::default(),
257        )
258        .unwrap();
259        assert_valid(&model);
260        assert!(model.sheet_face_names().is_empty(), "the plane is gone");
261        let chains = model.topology().wire_chains(wire).unwrap();
262        assert_eq!(chains.len(), 1);
263        assert_eq!(chains[0].len(), 4);
264        for &(edge, _) in &chains[0] {
265            let curve = &model.topology().get_edge(edge).unwrap().curve;
266            let mid = curve.evaluate(curve.mid_parameter().unwrap()).unwrap();
267            assert!(mid[2].could_be_equal(Scalar::from_f64(0.7)), "{mid:?}");
268            let r = mid[0].mul(mid[0]).add(mid[1].mul(mid[1])).sqrt().unwrap();
269            assert!((r.to_f64() - 0.5).abs() < 1e-9, "{r:?}");
270        }
271    }
272
273    /// A plane beside the body crosses nothing, and is refused.
274    #[test]
275    fn a_plane_beside_the_body_is_refused() {
276        let mut model = NamedModel::new();
277        let cube = block(&mut model, "a", [0.0; 3], [1.0; 3]);
278        let namer = Namer::new("section", "s").unwrap();
279        let plane = CoordinateSystem::world_at(v3(0.0, 0.0, 3.0));
280        let error = section_by_plane(
281            &mut model,
282            &namer,
283            Body::Solid(cube),
284            &plane,
285            RemeshParams::default(),
286        )
287        .unwrap_err();
288        assert!(error.root_message().contains("do not meet"), "{error:?}");
289        assert!(model.sheet_face_names().is_empty(), "the plane is gone");
290    }
291}
geop_lib_booleans/section.rs

geop_lib_booleans/
section.rs
