ogeom-offset 0.9.0

//! Pipes and lofts: the sweeps whose surfaces are already in the vocabulary.
//!
//! A circular profile along a straight spine is a cylinder; around a full
//! circle it is a torus; along an arc it is a torus segment, built from two
//! half-tube patches and two meridian caps so that every edge is a circle
//! with a closed-form chart. A ruled loft between two parallel sections is
//! walls of planes and cones: segment to segment gives the planar quad,
//! coaxial circle to circle gives the frustum the cone primitive already
//! builds, and a skew ruled wall between two segments is the bilinear patch
//! through its corners, exact. The sweeps that need surfaces with no closed
//! form (a free-form spine, a smoothed skin through many sections) are
//! fitted: the skin holds every section to a stated tolerance and refuses
//! when it cannot reach it.

use ogeom_algo::{
    Built, History, edge_vertices, make_cone, make_cylinder, make_edge, make_edge_between,
    make_face_with_pcurves, make_solid, make_torus, make_vertex, sew,
};
use ogeom_core::{OgeomResult, Tolerances, ogeom_bail};
use ogeom_geom::Curve3d as _;
use ogeom_geom::Transformable as _;
use ogeom_geom::{
    CircleCurve, Curve, Line2d, LineCurve, PlaneSurface, SurfaceGeometry, TorusSurface,
};
use ogeom_math::{Circle, Direction, Frame, Plane, Point, Point2, Torus, Transform, Vector};
use ogeom_topo::{EdgeData, EdgeRepr, Filter, Model, Shape, ShapeType, VertexData, explore};

/// Sweep a circular profile of `radius` along a spine edge.
///
/// A straight spine gives a cylinder, a full circular spine a torus, an arc
/// a torus segment. The history generates the solid from the spine.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the spine is
/// not straight or circular, the profile radius is not usable, or the tube
/// would swallow its own spine.
pub fn make_pipe(
    model: &mut Model,
    spine: &Shape,
    radius: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if !radius.is_finite() || radius <= tol.confusion() {
        ogeom_bail!(Construction, "a pipe of radius {radius} holds nothing");
    }
    let (curve, range) = {
        let Some(data) = model.node(spine).and_then(|n| n.data().as_edge()) else {
            ogeom_bail!(Construction, "a pipe runs along an edge");
        };
        let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
            ogeom_bail!(Construction, "the spine has no curve");
        };
        let Some(geometry) = model.geometry().curve(*curve) else {
            ogeom_bail!(Dangling, "curve is not in this model");
        };
        (geometry.clone(), *range)
    };
    let mut built = match &curve {
        Curve::Line(line) => {
            let start = curve.point_at(range.0, tol)?;
            let length = range.1 - range.0;
            let frame = Frame::about(start, line.axis().direction);
            make_cylinder(model, frame, radius, length, tol)?
        }
        Curve::Circle(c) => {
            let circle = c.circle();
            if radius >= circle.radius() - tol.confusion() {
                ogeom_bail!(
                    Construction,
                    "a tube of radius {radius} swallows its spine of radius {}",
                    circle.radius()
                );
            }
            let closed = curve
                .point_at(range.0, tol)?
                .distance(curve.point_at(range.1, tol)?)
                <= tol.confusion();
            if closed {
                make_torus(model, circle.frame(), circle.radius(), radius, tol)?
            } else {
                pipe_segment(model, circle, range, radius, tol)?
            }
        }
        _ => ogeom_bail!(
            Construction,
            "a pipe along a free-form spine needs the sweep-surface \
             machinery; see docs/PARITY.md, offset.sweeps"
        ),
    };
    built.history.generate(spine, built.shape.clone());
    Ok(built)
}

/// The torus segment: two half-tube patches and two meridian caps.
///
/// The tube circles are framed so their own parameter *is* the torus tube
/// angle, which makes every tube pcurve a vertical line in the chart and the
/// two patches the clean rectangles `v ∈ [0, π]` and `[π, 2π]`. The outer
/// equator is then the seam between the halves across the period (one edge,
/// two chart rows), which is exactly what [`ogeom_algo::attach_seam`] exists to
/// say.
fn pipe_segment(
    model: &mut Model,
    spine: Circle,
    range: (f64, f64),
    radius: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    let frame = spine.frame();
    let (x, y, z) = (frame.x().vector(), frame.y().vector(), frame.z().vector());
    let major = spine.radius();
    let radial = |u: f64| x * u.cos() + y * u.sin();
    let tangent = |u: f64| x * -u.sin() + y * u.cos();
    let tube_point = |u: f64, v: f64| {
        frame.origin() + radial(u) * radius.mul_add(v.cos(), major) + z * (radius * v.sin())
    };
    let pi = core::f64::consts::PI;
    let tau = core::f64::consts::TAU;

    let torus: SurfaceGeometry = TorusSurface::new(Torus::new(frame, major, radius, tol)?).into();
    let surface_id = model.geometry_mut().add_surface(torus);

    // Vertices at the tube's v = 0 and v = π points of each end.
    let ends = [range.0, range.1];
    let mut verts: Vec<Vec<Shape>> = Vec::new();
    for &u in &ends {
        verts.push(vec![
            make_vertex(model, tube_point(u, 0.0)).shape,
            make_vertex(model, tube_point(u, pi)).shape,
        ]);
    }

    // The tube circles at each end, split at v = 0 and v = π, framed so that
    // the circle's parameter equals the torus tube angle: `z` against the
    // spine tangent makes the frame's `y` the torus's own axis.
    let mut tube_arcs: Vec<Vec<Shape>> = Vec::new();
    for (k, &u) in ends.iter().enumerate() {
        let centre = frame.origin() + radial(u) * major;
        let circle = Circle::new(
            Frame::new(
                centre,
                Direction::new(-tangent(u), tol)?,
                Direction::new(radial(u), tol)?,
                tol,
            )?,
            radius,
            tol,
        )?;
        let curve = Curve::Circle(CircleCurve::new(circle));
        let arcs = vec![
            make_edge_between(
                model,
                curve.clone(),
                (0.0, pi),
                &verts[k][0],
                &verts[k][1],
                tol,
            )?
            .shape,
            make_edge_between(model, curve, (pi, tau), &verts[k][1], &verts[k][0], tol)?.shape,
        ];
        // In the chart both arcs run up the column at this end's angle.
        let column = Line2d::over(
            ogeom_math::Axis2::new(Point2::new(u, 0.0), ogeom_math::Direction2::Y),
            0.0,
            tau,
        )?;
        ogeom_algo::attach_pcurve(
            model,
            &arcs[0],
            column.into(),
            surface_id,
            ogeom_topo::Location::identity(),
            (0.0, pi),
        )?;
        ogeom_algo::attach_pcurve(
            model,
            &arcs[1],
            column.into(),
            surface_id,
            ogeom_topo::Location::identity(),
            (pi, tau),
        )?;
        tube_arcs.push(arcs);
    }

    // The long edges: the parallels at v = 0 and v = π, parameterized by the
    // spine's own angle.
    let parallel = |model: &mut Model, v: f64, from: &Shape, to: &Shape| -> OgeomResult<Shape> {
        let height = radius * v.sin();
        let ring = radius.mul_add(v.cos(), major);
        let circle = Circle::new(
            Frame::new(frame.origin() + z * height, frame.z(), frame.x(), tol)?,
            ring,
            tol,
        )?;
        let curve = Curve::Circle(CircleCurve::new(circle));
        Ok(make_edge_between(model, curve, range, from, to, tol)?.shape)
    };
    let row = |v: f64| -> OgeomResult<Line2d> {
        Line2d::over(
            ogeom_math::Axis2::new(Point2::new(0.0, v), ogeom_math::Direction2::X),
            range.0 - 1.0,
            range.1 + 1.0,
        )
    };
    let inner = parallel(model, pi, &verts[0][1], &verts[1][1])?;
    ogeom_algo::attach_pcurve(
        model,
        &inner,
        row(pi)?.into(),
        surface_id,
        ogeom_topo::Location::identity(),
        range,
    )?;
    // The outer equator bounds both halves across the period: v = 2π for its
    // forward use under the upper patch, v = 0 for its reversed use under the
    // lower: a seam, said as one.
    let outer = parallel(model, 0.0, &verts[0][0], &verts[1][0])?;
    ogeom_algo::attach_seam(
        model,
        &outer,
        row(tau)?.into(),
        row(0.0)?.into(),
        surface_id,
        ogeom_topo::Location::identity(),
        range,
    )?;

    // The two half-tube patches, on the one registered surface.
    let lower = {
        let wire = ogeom_algo::make_wire(
            model,
            &[
                tube_arcs[0][0].clone(),
                inner.clone(),
                tube_arcs[1][0].reversed(),
                outer.reversed(),
            ],
            tol,
        )?
        .shape;
        ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape
    };
    let upper = {
        let wire = ogeom_algo::make_wire(
            model,
            &[
                tube_arcs[0][1].clone(),
                outer.clone(),
                tube_arcs[1][1].reversed(),
                inner.reversed(),
            ],
            tol,
        )?
        .shape;
        ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape
    };

    // The meridian caps, their outward normals along the spine and away from
    // the material between the ends.
    let mut caps: Vec<Shape> = Vec::new();
    for (k, &u) in ends.iter().enumerate() {
        let outward = if k == 0 { -tangent(u) } else { tangent(u) };
        let centre = frame.origin() + radial(u) * major;
        let plane = Plane::through(centre, Direction::new(outward, tol)?);
        let reach = (major + radius) * 2.0;
        let surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        caps.push(
            make_face_with_pcurves(
                model,
                surface,
                &[vec![tube_arcs[k][0].clone(), tube_arcs[k][1].clone()]],
                tol,
            )?
            .shape,
        );
    }

    let faces = [lower, upper, caps[0].clone(), caps[1].clone()];
    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the pipe segment did not close");
    }
    make_solid(model, std::slice::from_ref(&sewn.shells[0]))
}

/// Loft two parallel closed sections into a solid, ruled.
///
/// Two coaxial circles give the cylinder or the cone frustum; two polygons
/// with the same corner count give planar walls. The sections pair edge by
/// edge in traversal order.
///
/// A wall between two segments that are not coplanar is the bilinear patch
/// through its four corners (the ruled surface between them, exact), so
/// sections may be turned against each other or differ in shape.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the sections
/// are not both circles or both polygons of the same count, or are not on
/// parallel planes.
pub fn make_loft(
    model: &mut Model,
    bottom: &Shape,
    top: &Shape,
    tol: Tolerances,
) -> OgeomResult<Built> {
    // A vertex for either section is the loft to a point: a cone over a
    // circle, an exact pyramid over a polygon.
    match (model.kind_of(bottom)?, model.kind_of(top)?) {
        (ShapeType::Wire, ShapeType::Vertex) => {
            return loft_to_point(model, bottom, top, tol);
        }
        (ShapeType::Vertex, ShapeType::Wire) => {
            let mut built = loft_to_point(model, top, bottom, tol)?;
            // The apex was named first: the same solid, the history the same.
            built.history.generate(bottom, built.shape.clone());
            return Ok(built);
        }
        _ => {}
    }
    for wire in [bottom, top] {
        if model.kind_of(wire)? != ShapeType::Wire {
            ogeom_bail!(Construction, "a loft runs between wires");
        }
        if !ogeom_algo::is_wire_closed(model, wire, tol)? {
            ogeom_bail!(Construction, "a loft section must be closed");
        }
    }
    let circle_of = |model: &Model, wire: &Shape| -> OgeomResult<Option<Circle>> {
        let edges = explore(model, wire, Filter::OfType(ShapeType::Edge))?;
        if edges.len() != 1 {
            return Ok(None);
        }
        let Some(data) = model.node(&edges[0]).and_then(|n| n.data().as_edge()) else {
            return Ok(None);
        };
        let Some(EdgeRepr::Curve3d { curve, .. }) = data.curve3d() else {
            return Ok(None);
        };
        match model.geometry().curve(*curve) {
            Some(Curve::Circle(c)) => Ok(Some(c.circle())),
            _ => Ok(None),
        }
    };

    if let (Some(lower), Some(upper)) = (circle_of(model, bottom)?, circle_of(model, top)?) {
        // Coaxial circles: the revolved primitives already build these.
        let axis = lower.frame().z().vector();
        let rise = upper.centre() - lower.centre();
        let height = rise.dot(axis);
        if rise.cross(axis).magnitude() > tol.confusion() * 10.0 || height.abs() <= tol.confusion()
        {
            ogeom_bail!(
                Construction,
                "lofted circles must be coaxial on parallel planes; the \
                 oblique loft needs the sweep machinery; see the deferred \
                 table"
            );
        }
        let frame = if height > 0.0 {
            Frame::new(lower.centre(), lower.frame().z(), lower.frame().x(), tol)?
        } else {
            Frame::new(
                lower.centre(),
                lower.frame().z().reversed(),
                lower.frame().x(),
                tol,
            )?
        };
        let mut built = if (lower.radius() - upper.radius()).abs() <= tol.confusion() {
            make_cylinder(model, frame, lower.radius(), height.abs(), tol)?
        } else {
            make_cone(
                model,
                frame,
                lower.radius(),
                upper.radius(),
                height.abs(),
                tol,
            )?
        };
        built.history.generate(bottom, built.shape.clone());
        built.history.generate(top, built.shape.clone());
        return Ok(built);
    }

    // Polygonal sections: matched corners, planar walls.
    let corners_of = |model: &Model, wire: &Shape| -> OgeomResult<Vec<Point>> {
        let mut out = Vec::new();
        for edge in explore(model, wire, Filter::OfType(ShapeType::Edge))? {
            let Some(data) = model.node(&edge).and_then(|n| n.data().as_edge()) else {
                ogeom_bail!(Construction, "a section edge holds no data");
            };
            let Some(EdgeRepr::Curve3d { curve, .. }) = data.curve3d() else {
                ogeom_bail!(Construction, "a section edge has no curve");
            };
            let Some(Curve::Line(_)) = model.geometry().curve(*curve) else {
                ogeom_bail!(
                    Construction,
                    "a mixed or curved section needs the skinning machinery; see \
                     docs/PARITY.md, offset.sweeps"
                );
            };
            let Some((sv, _)) = edge_vertices(model, &edge)? else {
                ogeom_bail!(Construction, "a section edge has no vertices");
            };
            let Some(data) = model.node(&sv).and_then(|n| n.data().as_vertex()) else {
                ogeom_bail!(Construction, "a section vertex holds no point");
            };
            out.push(sv.transform(model.datums())?.apply(data.point));
        }
        Ok(out)
    };
    let low = corners_of(model, bottom)?;
    let high = corners_of(model, top)?;
    if low.len() != high.len() {
        ogeom_bail!(
            Construction,
            "lofted sections must have the same corner count, found {} and {}",
            low.len(),
            high.len()
        );
    }
    let n = low.len();
    let centroid = {
        let mut c = Vector::new(0.0, 0.0, 0.0);
        for p in low.iter().chain(high.iter()) {
            c += p.to_vector();
        }
        #[allow(clippy::cast_precision_loss)]
        let count = 2.0 * n as f64;
        Point::from_vector(c / count)
    };

    // Shared vertices and edges, then walls and caps referencing them.
    let vl: Vec<Shape> = low.iter().map(|p| make_vertex(model, *p).shape).collect();
    let vh: Vec<Shape> = high.iter().map(|p| make_vertex(model, *p).shape).collect();
    let seg = |model: &mut Model, a: (&Shape, Point), b: (&Shape, Point)| -> OgeomResult<Shape> {
        let line = LineCurve::segment(a.1, b.1, tol)?;
        let curve = Curve::Line(line);
        let domain = curve.domain();
        Ok(make_edge_between(model, curve, domain, a.0, b.0, tol)?.shape)
    };
    let mut low_edges = Vec::with_capacity(n);
    let mut high_edges = Vec::with_capacity(n);
    let mut rails = Vec::with_capacity(n);
    for i in 0..n {
        let j = (i + 1) % n;
        low_edges.push(seg(model, (&vl[i], low[i]), (&vl[j], low[j]))?);
        high_edges.push(seg(model, (&vh[i], high[i]), (&vh[j], high[j]))?);
        rails.push(seg(model, (&vl[i], low[i]), (&vh[i], high[i]))?);
    }

    let planar = |model: &mut Model, corners: &[Point], edges: Vec<Shape>| -> OgeomResult<Shape> {
        let normal = {
            let mut n = (corners[1] - corners[0]).cross(corners[2] - corners[0]);
            let m = n.magnitude();
            if m <= tol.confusion() {
                ogeom_bail!(Construction, "a loft wall is degenerate");
            }
            n /= m;
            if n.dot(corners[0] - centroid) < 0.0 {
                -n
            } else {
                n
            }
        };
        let skew = corners.iter().any(|p| {
            Plane::through(
                corners[0],
                Direction::new(normal, tol).unwrap_or(Direction::Z),
            )
            .distance_to(*p)
                > tol.confusion() * 10.0
        });
        if skew {
            // A wall between two segments that do not lie in one plane is
            // the ruled surface between them, and that is exact: bilinear
            // in the four corners, a B-spline of degree one each way.
            if corners.len() != 4 || edges.len() != 4 {
                ogeom_bail!(Construction, "a skew ruled wall has four corners");
            }
            return bilinear_wall(model, corners, &edges, centroid, tol);
        }
        let plane = Plane::through(corners[0], Direction::new(normal, tol)?);
        let mut reach = 1.0_f64;
        for p in corners {
            reach = reach.max(p.distance(corners[0]) * 2.0);
        }
        let surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        Ok(make_face_with_pcurves(model, surface, &[edges], tol)?.shape)
    };

    let mut faces: Vec<Shape> = Vec::with_capacity(n + 2);
    for i in 0..n {
        let j = (i + 1) % n;
        faces.push(planar(
            model,
            &[low[i], low[j], high[j], high[i]],
            vec![
                low_edges[i].clone(),
                rails[j].clone(),
                high_edges[i].reversed(),
                rails[i].reversed(),
            ],
        )?);
    }
    faces.push(planar(model, &low, low_edges.clone())?);
    faces.push(planar(model, &high, high_edges.clone())?);

    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the loft did not close");
    }
    let mut built = make_solid(model, std::slice::from_ref(&sewn.shells[0]))?;
    built.history.generate(bottom, built.shape.clone());
    built.history.generate(top, built.shape.clone());
    Ok(built)
}

/// The ruled wall between two straight segments that are not coplanar: the
/// bilinear patch through its four corners, exact, as a B-spline of degree
/// one each way.
///
/// `corners` run round the wall (low start, low end, high end, high
/// start) and `edges` walk them in that order, the third and fourth
/// reversed as the caller's wire has them. Each edge's pcurve is the chart
/// side it lies on, parameterized by the edge's own range so the two agree
/// point for point; the wall faces away from `centroid`.
fn bilinear_wall(
    model: &mut Model,
    corners: &[Point],
    edges: &[Shape],
    centroid: Point,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    use ogeom_geom::Surface as _;
    let (p00, p10, p11, p01) = (corners[0], corners[1], corners[2], corners[3]);
    let grid = ogeom_math::ControlGrid::new(vec![p00, p01, p10, p11], 2, 2)?;
    let line = ogeom_math::KnotVector::clamped_uniform(1, 2)?;
    let patch = ogeom_geom::BSplineSurface::new(line.clone(), line, &grid, tol)?;
    let outward = {
        let (du, dv) = patch.d1_at(0.5, 0.5, tol)?;
        let centre = patch.point_at(0.5, 0.5, tol)?;
        du.cross(dv).dot(centre - centroid) >= 0.0
    };
    let surface_id = model
        .geometry_mut()
        .add_surface(SurfaceGeometry::BSpline(patch));

    // Each edge's chart side, run in the edge's own direction over the
    // edge's own parameter range: a degree-one B-spline in the chart,
    // which is what makes the pcurve the edge's equal parameter for
    // parameter whatever length the edge has.
    let sides: [(Point2, Point2); 4] = [
        (Point2::new(0.0, 0.0), Point2::new(1.0, 0.0)),
        (Point2::new(1.0, 0.0), Point2::new(1.0, 1.0)),
        (Point2::new(0.0, 1.0), Point2::new(1.0, 1.0)),
        (Point2::new(0.0, 0.0), Point2::new(0.0, 1.0)),
    ];
    for (edge, (from, to)) in edges.iter().zip(sides) {
        let range = {
            let Some(node) = model.node(edge) else {
                ogeom_bail!(Dangling, "a loft edge is not in this model");
            };
            let Some(data) = node.data().as_edge() else {
                ogeom_bail!(Construction, "a loft edge holds no edge data");
            };
            let Some(EdgeRepr::Curve3d { range, .. }) = data.curve3d() else {
                ogeom_bail!(Construction, "a loft edge has no curve");
            };
            *range
        };
        let knots = ogeom_math::KnotVector::new(vec![range.0, range.0, range.1, range.1], 1)?;
        let pcurve = ogeom_geom::BSpline2d::new(knots, vec![from, to], tol)?;
        ogeom_algo::attach_pcurve(
            model,
            edge,
            pcurve.into(),
            surface_id,
            ogeom_topo::Location::identity(),
            range,
        )?;
    }
    let wire = ogeom_algo::make_wire(model, edges, tol)?.shape;
    let face = ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape;
    Ok(if outward { face } else { face.reversed() })
}

/// A skinned wall and the pieces a caller needs to close it: the rings at
/// both ends, their exact border curves off the control net, and the chart's
/// `u` window the ring pcurves span.
struct SkinnedWall {
    face: Shape,
    ring0: Shape,
    ring1: Shape,
    curve0: ogeom_geom::Curve,
    curve1: ogeom_geom::Curve,
    u_dom: (f64, f64),
}

/// The wall of a skin over a grid of section samples, closed the way round.
///
/// The wall is [`ogeom_geom::fit::fit_surface_grid`]'s surface with each row's
/// first sample repeated at its end: the row fits pin their ends, so the two
/// border control columns are *equal* and the seam closes exactly, not
/// within tolerance. The border iso-curves come straight off the control
/// net (the v-borders are the fitted sections, planar whenever the
/// sections are, which is what lets the caps be planes), and every pcurve
/// is an iso line in the fitted chart, same-parameter by construction.
fn skinned_wall(
    model: &mut Model,
    rows: &[Vec<Point>],
    shared: (Option<&Shape>, Option<&Shape>),
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<SkinnedWall> {
    use ogeom_geom::Surface as _;
    let mut closed_rows: Vec<Vec<Point>> = Vec::with_capacity(rows.len());
    for row in rows {
        let mut r = row.clone();
        r.push(row[0]);
        closed_rows.push(r);
    }
    let fitted = ogeom_geom::fit::fit_surface_grid(&closed_rows, 3, tolerance, tol)?;
    if !fitted.met {
        ogeom_bail!(
            NotDone,
            "the skin reached {} against a target of {tolerance}",
            fitted.error
        );
    }
    let surface = fitted.curve;
    let (u_knots, v_knots) = (surface.u_knots().clone(), surface.v_knots().clone());
    let (k, l, net) = {
        let grid = surface.grid();
        let net: Vec<Point> = grid.points().iter().map(|w| (*w).point()).collect();
        (grid.u_count(), grid.v_count(), net)
    };
    let point_at = |i: usize, j: usize| -> Point { net[i * l + j] };
    let (u_dom, v_dom) = surface.domain();

    // Border curves straight off the net: v-borders are the end sections,
    // the u-border is the seam.
    let border_v = |j: usize| -> OgeomResult<ogeom_geom::Curve> {
        let control: Vec<Point> = (0..k).map(|i| point_at(i, j)).collect();
        Ok(ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            u_knots.clone(),
            control,
            tol,
        )?))
    };
    let seam_curve = {
        let control: Vec<Point> = (0..l).map(|j| point_at(0, j)).collect();
        ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            v_knots.clone(),
            control,
            tol,
        )?)
    };

    let surface_geo: SurfaceGeometry = surface.into();
    let surface_id = model.geometry_mut().add_surface(surface_geo.clone());

    // An end ring the neighbouring wall already built is adopted, not
    // refitted (see `adopt_border`).
    let slack = fitted.error + tol.confusion();
    let ring_of = |model: &mut Model,
                   j: usize,
                   given: Option<&Shape>|
     -> OgeomResult<(Shape, ogeom_geom::Curve)> {
        match given {
            Some(edge) => {
                adopt_border(model, edge, &surface_geo, slack, tol)?;
                Ok((edge.clone(), spine_curve_of(model, edge)?.0))
            }
            None => {
                let curve = border_v(j)?;
                Ok((make_edge(model, curve.clone(), u_dom, tol)?.shape, curve))
            }
        }
    };
    let (ring0, curve0) = ring_of(model, 0, shared.0)?;
    let (ring1, curve1) = ring_of(model, l - 1, shared.1)?;
    let anchor0 = ogeom_algo::edge_vertices(model, &ring0)?
        .map(|(a, _)| a)
        .ok_or_else(|| ogeom_core::ogeom_err!(Construction, "a skinned ring has no vertex"))?;
    let anchor1 = ogeom_algo::edge_vertices(model, &ring1)?
        .map(|(a, _)| a)
        .ok_or_else(|| ogeom_core::ogeom_err!(Construction, "a skinned ring has no vertex"))?;
    let seam = make_edge_between(model, seam_curve, v_dom, &anchor0, &anchor1, tol)?.shape;

    // Pcurves: rows for the rings, both columns for the seam.
    let row_line = |v: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(0.0, v), ogeom_math::Direction2::X),
            u_dom.0 - 1.0,
            u_dom.1 + 1.0,
        )?
        .into())
    };
    let column_line = |u: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(u, 0.0), ogeom_math::Direction2::Y),
            v_dom.0 - 1.0,
            v_dom.1 + 1.0,
        )?
        .into())
    };
    ogeom_algo::attach_pcurve(
        model,
        &ring0,
        row_line(v_dom.0)?,
        surface_id,
        ogeom_topo::Location::identity(),
        u_dom,
    )?;
    ogeom_algo::attach_pcurve(
        model,
        &ring1,
        row_line(v_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        u_dom,
    )?;
    ogeom_algo::attach_seam(
        model,
        &seam,
        column_line(u_dom.0)?,
        column_line(u_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        v_dom,
    )?;

    let wall = {
        let wire = ogeom_algo::make_wire(
            model,
            &[
                ring0.clone(),
                seam.clone(),
                ring1.reversed(),
                seam.reversed(),
            ],
            tol,
        )?
        .shape;
        let face =
            ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape;
        // Outward by measurement at the middle of the skin: away from the
        // centre of the section the point lies on. The centre of the whole
        // skin will not do: along a bent spine it lies inside the bend,
        // on the far side of the inner wall.
        let mid_u = f64::midpoint(u_dom.0, u_dom.1);
        let mid_v = f64::midpoint(v_dom.0, v_dom.1);
        let s_mid = surface_geo.point_at(mid_u, mid_v, tol)?;
        let (du, dv) = surface_geo.d1_at(mid_u, mid_v, tol)?;
        if du.cross(dv).dot(s_mid - section_centre(rows, s_mid)) >= 0.0 {
            face
        } else {
            face.reversed()
        }
    };
    Ok(SkinnedWall {
        face: wall,
        ring0,
        ring1,
        curve0,
        curve1,
        u_dom,
    })
}

/// The centre of the section (row) passing nearest `at`: the mean of its
/// points.
fn section_centre(rows: &[Vec<Point>], at: Point) -> Point {
    let nearest = rows.iter().filter(|row| !row.is_empty()).min_by(|a, b| {
        let d = |row: &Vec<Point>| {
            row.iter()
                .map(|p| p.distance(at))
                .fold(f64::INFINITY, f64::min)
        };
        d(a).total_cmp(&d(b))
    });
    let Some(row) = nearest else {
        return at;
    };
    let sum = row
        .iter()
        .fold(Vector::new(0.0, 0.0, 0.0), |acc, p| acc + p.to_vector());
    #[allow(clippy::cast_precision_loss)]
    Point::from_vector(sum / row.len() as f64)
}

/// A strip closed the *long* way: open across its own width, a smooth loop
/// along the sweep: one face of a faceted ring, [`skinned_wall`]'s
/// construction with the chart's roles swapped and the loop made C1 by
/// [`ogeom_geom::fit::fit_surface_grid_closed_v`]. The rails are the two
/// closed border loops; the seam is one station's column, used twice.
fn skinned_ring_strip(
    model: &mut Model,
    rows: &[Vec<Point>],
    outward_hint: Point,
    shared: [Option<&Shape>; 2],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<(Shape, Shape, Shape)> {
    use ogeom_geom::Surface as _;
    // The loop: first row repeated at the end, as the closed fit demands.
    let mut looped: Vec<Vec<Point>> = rows.to_vec();
    looped.push(rows[0].clone());
    let fitted = ogeom_geom::fit::fit_surface_grid_closed_v(&looped, 3, tolerance, tol)?;
    if !fitted.met {
        ogeom_bail!(
            NotDone,
            "the ring strip reached {} against a target of {tolerance}",
            fitted.error
        );
    }
    let surface = fitted.curve;
    let (u_knots, v_knots) = (surface.u_knots().clone(), surface.v_knots().clone());
    let (k, l, net) = {
        let grid = surface.grid();
        let net: Vec<Point> = grid.points().iter().map(|w| (*w).point()).collect();
        (grid.u_count(), grid.v_count(), net)
    };
    let point_at = |i: usize, j: usize| -> Point { net[i * l + j] };
    let (u_dom, v_dom) = surface.domain();

    // The chart's roles, straight: `u` runs across the strip (open), `v`
    // around the loop (closed). The rails are v-curves (the closed border
    // loops at the two u-borders), and the seam is the u-row at the loop's
    // join, bounding the chart twice as every seam does.
    let rail_curve = |i: usize| -> OgeomResult<ogeom_geom::Curve> {
        let control: Vec<Point> = (0..l).map(|j| point_at(i, j)).collect();
        Ok(ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            v_knots.clone(),
            control,
            tol,
        )?))
    };
    let seam_curve = {
        let control: Vec<Point> = (0..k).map(|i| point_at(i, 0)).collect();
        ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            u_knots.clone(),
            control,
            tol,
        )?)
    };
    let surface_geo: SurfaceGeometry = surface.into();
    let surface_id = model.geometry_mut().add_surface(surface_geo.clone());

    // A corner loop shared with the neighbouring strip is one edge for
    // both: the neighbour built it from its own fit, and this strip's
    // border is another fit of the same loop, so the edge widens to how
    // far it honestly sits from this surface. Two independent fits of one
    // loop can disagree by more than either fit's own error (the sew refuses
    // that gap under a frame that turns fast), and one edge cannot disagree
    // with itself.
    let slack = fitted.error + tol.confusion();
    let rail_of = |model: &mut Model, i: usize, given: Option<&Shape>| -> OgeomResult<Shape> {
        let Some(edge) = given else {
            let edge = make_edge(model, rail_curve(i)?, v_dom, tol)?.shape;
            model.widen(&edge, ogeom_core::Tolerance::new(slack)?)?;
            return Ok(edge);
        };
        let (curve, range) = spine_curve_of(model, edge)?;
        let mut off: f64 = 0.0;
        for step in 0..=32 {
            #[allow(clippy::cast_precision_loss)]
            let t = range.0 + (range.1 - range.0) * (step as f64) / 32.0;
            let p = curve.point_at(t, tol)?;
            off = off.max(ogeom_algo::project_on_surface(&surface_geo, p, 16, tol)?.distance);
        }
        model.widen(edge, ogeom_core::Tolerance::new(off + slack)?)?;
        if let Some((a, b)) = ogeom_algo::edge_vertices(model, edge)? {
            for v in [&a, &b] {
                model.widen(v, ogeom_core::Tolerance::new(off + slack)?)?;
            }
        }
        Ok(edge.clone())
    };
    let rail0 = rail_of(model, 0, shared[0])?;
    let rail1 = rail_of(model, k - 1, shared[1])?;
    let anchor0 = ogeom_algo::edge_vertices(model, &rail0)?
        .map(|(a, _)| a)
        .ok_or_else(|| ogeom_core::ogeom_err!(Construction, "a strip rail has no vertex"))?;
    let anchor1 = ogeom_algo::edge_vertices(model, &rail1)?
        .map(|(a, _)| a)
        .ok_or_else(|| ogeom_core::ogeom_err!(Construction, "a strip rail has no vertex"))?;
    let seam = make_edge_between(model, seam_curve, u_dom, &anchor0, &anchor1, tol)?.shape;

    let row_line = |v: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(0.0, v), ogeom_math::Direction2::X),
            u_dom.0 - 1.0,
            u_dom.1 + 1.0,
        )?
        .into())
    };
    let column_line = |u: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(u, 0.0), ogeom_math::Direction2::Y),
            v_dom.0 - 1.0,
            v_dom.1 + 1.0,
        )?
        .into())
    };
    ogeom_algo::attach_pcurve(
        model,
        &rail0,
        column_line(u_dom.0)?,
        surface_id,
        ogeom_topo::Location::identity(),
        v_dom,
    )?;
    ogeom_algo::attach_pcurve(
        model,
        &rail1,
        column_line(u_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        v_dom,
    )?;
    ogeom_algo::attach_seam(
        model,
        &seam,
        row_line(v_dom.0)?,
        row_line(v_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        u_dom,
    )?;
    let wire = ogeom_algo::make_wire(
        model,
        &[
            rail0.clone(),
            seam.clone(),
            rail1.reversed(),
            seam.reversed(),
        ],
        tol,
    )?
    .shape;
    let face = ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape;
    let mid_u = f64::midpoint(u_dom.0, u_dom.1);
    let mid_v = f64::midpoint(v_dom.0, v_dom.1);
    let s_mid = surface_geo.point_at(mid_u, mid_v, tol)?;
    let (du, dv) = surface_geo.d1_at(mid_u, mid_v, tol)?;
    let face = if du.cross(dv).dot(s_mid - outward_hint) >= 0.0 {
        face
    } else {
        face.reversed()
    };
    Ok((face, rail0, rail1))
}

/// A border edge the neighbouring skin built from its own fit, adopted by
/// this skin: the edge widens to how far it honestly sits from `surface`
/// plus this fit's own slack, its vertices with it. Two independent fits of
/// one row can disagree by more than either fit's own error, and one edge
/// cannot disagree with itself.
fn adopt_border(
    model: &mut Model,
    edge: &Shape,
    surface: &SurfaceGeometry,
    slack: f64,
    tol: Tolerances,
) -> OgeomResult<()> {
    let (curve, range) = spine_curve_of(model, edge)?;
    let mut off: f64 = 0.0;
    for step in 0..=32 {
        #[allow(clippy::cast_precision_loss)]
        let t = range.0 + (range.1 - range.0) * (step as f64) / 32.0;
        let p = curve.point_at(t, tol)?;
        off = off.max(ogeom_algo::project_on_surface(surface, p, 16, tol)?.distance);
    }
    model.widen(edge, ogeom_core::Tolerance::new(off + slack)?)?;
    if let Some((a, b)) = ogeom_algo::edge_vertices(model, edge)? {
        for v in [&a, &b] {
            model.widen(v, ogeom_core::Tolerance::new(off + slack)?)?;
        }
    }
    Ok(())
}

/// An adopted border's image on a surface it was not fitted on: the
/// border's own points, each read off the surface at its nearest point,
/// fitted at the border's own parameters so the image is same-parameter
/// with it. Returned with the parameter range it spans.
fn adopted_image(
    model: &Model,
    edge: &Shape,
    surface: &SurfaceGeometry,
    tol: Tolerances,
) -> OgeomResult<(ogeom_geom::PlanarCurve, (f64, f64))> {
    const SAMPLES: u32 = 64;
    let (curve, range) = spine_curve_of(model, edge)?;
    let mut params = Vec::with_capacity(SAMPLES as usize + 1);
    let mut image = Vec::with_capacity(SAMPLES as usize + 1);
    let mut guess: Option<(f64, f64)> = None;
    for step in 0..=SAMPLES {
        let t = range.0 + (range.1 - range.0) * f64::from(step) / f64::from(SAMPLES);
        let p = curve.point_at(t, tol)?;
        let foot = match guess {
            Some(g) => ogeom_algo::project_on_surface_from(surface, p, g, tol)?,
            None => ogeom_algo::project_on_surface(surface, p, 16, tol)?,
        };
        guess = Some(foot.parameters);
        params.push(t);
        image.push(Point2::new(foot.parameters.0, foot.parameters.1));
    }
    let fitted = ogeom_geom::fit::fit_points_2d_at(&params, &image, 3, tol.confusion(), tol)?;
    Ok((fitted.curve.into(), range))
}

/// A solid skinned over a grid of section samples: [`skinned_wall`] with a
/// planar cap over each end ring.
/// How a skinned solid's end is closed.
#[derive(Debug, Clone, Copy)]
enum EndCap {
    /// The section is planar: a plane face, its normal pointing out.
    Plane(Vector),
    /// The section is not: a patch skinned from the ring down to a point
    /// inside it, sharing the wall's ring edge.
    Skinned,
}

fn skinned_solid(
    model: &mut Model,
    rows: &[Vec<Point>],
    caps: (EndCap, EndCap),
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    let wall = skinned_wall(model, rows, (None, None), tolerance, tol)?;
    let u_dom = wall.u_dom;

    let cap = |model: &mut Model,
               ring: &Shape,
               curve: ogeom_geom::Curve,
               outward: Vector|
     -> OgeomResult<Shape> {
        let at = curve.point_at(u_dom.0, tol)?;
        let plane = Plane::through(at, Direction::new(outward, tol)?);
        let mut reach = 1.0_f64;
        for t in 0..8 {
            let p = curve.point_at(u_dom.0 + (u_dom.1 - u_dom.0) * f64::from(t) / 8.0, tol)?;
            reach = reach.max(p.distance(at) * 2.0);
        }
        let cap_surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        let wire = ogeom_algo::make_wire(model, std::slice::from_ref(ring), tol)?.shape;
        let face =
            ogeom_algo::make_face(model, cap_surface.clone(), std::slice::from_ref(&wire), tol)?
                .shape;
        let id = {
            let Some(node) = model.node(&face) else {
                ogeom_bail!(Dangling, "the cap just built is not in this model");
            };
            let ogeom_topo::NodeData::Face(data) = node.data() else {
                ogeom_bail!(Construction, "the cap holds no face data");
            };
            data.surface
        };
        let Some(pcurve) = ogeom_intersect::exact_pcurve_of(&curve, &cap_surface, tol) else {
            ogeom_bail!(Construction, "a cap edge has no closed-form pcurve");
        };
        ogeom_algo::attach_pcurve(
            model,
            ring,
            pcurve,
            id,
            ogeom_topo::Location::identity(),
            u_dom,
        )?;
        Ok(face)
    };
    let close =
        |model: &mut Model, end: EndCap, ring: &Shape, curve: &ogeom_geom::Curve, row: &[Point]| {
            match end {
                EndCap::Plane(outward) => cap(model, ring, curve.clone(), outward),
                EndCap::Skinned => {
                    // The ring, a row halfway in, and the point the rest of
                    // the ring's rows collapse to: the section's own centroid,
                    // which a closed section winds round.
                    let apex = centroid_of(std::slice::from_ref(&row.to_vec()));
                    let half: Vec<Point> = row
                        .iter()
                        .map(|p| Point::from_vector((p.to_vector() + apex.to_vector()) * 0.5))
                        .collect();
                    let rows = [row.to_vec(), half, vec![apex; row.len()]];
                    Ok(apex_patch(model, &rows, Some(ring), tolerance, tol)?.0)
                }
            }
        };
    let cap0 = close(model, caps.0, &wall.ring0, &wall.curve0, &rows[0])?;
    let cap1 = close(
        model,
        caps.1,
        &wall.ring1,
        &wall.curve1,
        &rows[rows.len() - 1],
    )?;

    let faces = [wall.face, cap0, cap1];
    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the skinned solid did not close");
    }
    make_solid(model, std::slice::from_ref(&sewn.shells[0]))
}

/// A patch skinned from a ring down to a point: [`skinned_wall`]'s
/// construction with the top ring replaced by the apex: a degenerate
/// edge on one vertex, bounding the chart's whole top row the way a cone's
/// apex bounds a countersink. The ring edge is adopted from `shared`
/// where a neighbour already built it, and the face is turned to point
/// away from the section it passes through. Returns the face and its ring
/// edge.
fn apex_patch(
    model: &mut Model,
    rows: &[Vec<Point>],
    shared: Option<&Shape>,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<(Shape, Shape)> {
    use ogeom_geom::Surface as _;
    let mut closed_rows: Vec<Vec<Point>> = Vec::with_capacity(rows.len());
    for row in rows {
        let mut r = row.clone();
        r.push(row[0]);
        closed_rows.push(r);
    }
    let fitted = ogeom_geom::fit::fit_surface_grid(&closed_rows, 3, tolerance, tol)?;
    if !fitted.met {
        ogeom_bail!(
            NotDone,
            "the skin reached {} against a target of {tolerance}",
            fitted.error
        );
    }
    let surface = fitted.curve;
    let (u_knots, v_knots) = (surface.u_knots().clone(), surface.v_knots().clone());
    let (k, l, net) = {
        let grid = surface.grid();
        let net: Vec<Point> = grid.points().iter().map(|w| (*w).point()).collect();
        (grid.u_count(), grid.v_count(), net)
    };
    let point_at = |i: usize, j: usize| -> Point { net[i * l + j] };
    let (u_dom, v_dom) = surface.domain();
    let apex = rows[rows.len() - 1][0];

    let ring_curve = {
        let control: Vec<Point> = (0..k).map(|i| point_at(i, 0)).collect();
        ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            u_knots.clone(),
            control,
            tol,
        )?)
    };
    let seam_curve = {
        let control: Vec<Point> = (0..l).map(|j| point_at(0, j)).collect();
        ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            v_knots.clone(),
            control,
            tol,
        )?)
    };
    let surface_geo: SurfaceGeometry = surface.into();
    let surface_id = model.geometry_mut().add_surface(surface_geo.clone());

    // The ring a neighbour built is adopted, not refitted (see `adopt_border`).
    let ring0 = match shared {
        Some(edge) => {
            adopt_border(
                model,
                edge,
                &surface_geo,
                fitted.error + tol.confusion(),
                tol,
            )?;
            edge.clone()
        }
        None => make_edge(model, ring_curve.clone(), u_dom, tol)?.shape,
    };
    let anchor0 = ogeom_algo::edge_vertices(model, &ring0)?
        .map(|(a, _)| a)
        .ok_or_else(|| ogeom_core::ogeom_err!(Construction, "a skinned ring has no vertex"))?;
    let apex_vertex = model.add_vertex(VertexData::new(apex));
    let apex_edge = {
        let mut data = EdgeData::new();
        data.degenerate = true;
        model.add_edge(data, &[apex_vertex.clone(), apex_vertex.clone()])?
    };
    let seam = make_edge_between(model, seam_curve, v_dom, &anchor0, &apex_vertex, tol)?.shape;

    let row_line = |v: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(0.0, v), ogeom_math::Direction2::X),
            u_dom.0 - 1.0,
            u_dom.1 + 1.0,
        )?
        .into())
    };
    let column_line = |u: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(u, 0.0), ogeom_math::Direction2::Y),
            v_dom.0 - 1.0,
            v_dom.1 + 1.0,
        )?
        .into())
    };
    ogeom_algo::attach_pcurve(
        model,
        &ring0,
        row_line(v_dom.0)?,
        surface_id,
        ogeom_topo::Location::identity(),
        u_dom,
    )?;
    // The apex bounds the chart's whole top row while covering no distance:
    // the degenerate edge carries the row's pcurve, exactly as a cone's apex
    // does after the reader synthesises it.
    ogeom_algo::attach_pcurve(
        model,
        &apex_edge,
        row_line(v_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        u_dom,
    )?;
    ogeom_algo::attach_seam(
        model,
        &seam,
        column_line(u_dom.0)?,
        column_line(u_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        v_dom,
    )?;

    let wire = ogeom_algo::make_wire(
        model,
        &[
            ring0.clone(),
            seam.clone(),
            apex_edge.reversed(),
            seam.reversed(),
        ],
        tol,
    )?
    .shape;
    let face = ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape;
    let mid_u = f64::midpoint(u_dom.0, u_dom.1);
    let mid_v = f64::midpoint(v_dom.0, v_dom.1);
    let s_mid = surface_geo.point_at(mid_u, mid_v, tol)?;
    let (du, dv) = surface_geo.d1_at(mid_u, mid_v, tol)?;
    // Inside is judged from the section the patch passes through at that
    // row, not from the whole skin's centroid: a skin that bends puts that
    // centroid outside itself, in the crook of the bend.
    let local = {
        let mut sum = Vector::ZERO;
        for k in 0..16 {
            let u = u_dom.0 + (u_dom.1 - u_dom.0) * f64::from(k) / 16.0;
            sum += surface_geo.point_at(u, mid_v, tol)?.to_vector();
        }
        Point::from_vector(sum / 16.0)
    };
    let face = if du.cross(dv).dot(s_mid - local) >= 0.0 {
        face
    } else {
        face.reversed()
    };
    Ok((face, ring0))
}

/// The mean of every point in every row.
fn centroid_of(rows: &[Vec<Point>]) -> Point {
    let mut c = Vector::new(0.0, 0.0, 0.0);
    let mut n = 0.0;
    for row in rows {
        for p in row {
            c += p.to_vector();
            n += 1.0;
        }
    }
    Point::from_vector(c / n)
}

/// A solid skinned down to a point: [`apex_patch`] for the wall, and one
/// cap at the open end; the apex end closes by construction.
fn skinned_solid_to_apex(
    model: &mut Model,
    rows: &[Vec<Point>],
    cap_outward: Vector,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    use ogeom_geom::Curve3d as _;
    let (wall, ring0) = apex_patch(model, rows, None, tolerance, tol)?;
    let (ring_curve, u_dom) = {
        let (curve, range) = spine_curve_of(model, &ring0)?;
        (curve, range)
    };

    // One cap, on the open end; the machinery is skinned_solid's, inlined
    // for the single ring.
    let cap = {
        let at = ring_curve.point_at(u_dom.0, tol)?;
        let plane = Plane::through(at, Direction::new(cap_outward, tol)?);
        let mut reach = 1.0_f64;
        for t in 0..8 {
            let p = ring_curve.point_at(u_dom.0 + (u_dom.1 - u_dom.0) * f64::from(t) / 8.0, tol)?;
            reach = reach.max(p.distance(at) * 2.0);
        }
        let cap_surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        let wire = ogeom_algo::make_wire(model, std::slice::from_ref(&ring0), tol)?.shape;
        let face =
            ogeom_algo::make_face(model, cap_surface.clone(), std::slice::from_ref(&wire), tol)?
                .shape;
        let id = {
            let Some(node) = model.node(&face) else {
                ogeom_bail!(Dangling, "the cap just built is not in this model");
            };
            let ogeom_topo::NodeData::Face(data) = node.data() else {
                ogeom_bail!(Construction, "the cap holds no face data");
            };
            data.surface
        };
        let Some(pcurve) = ogeom_intersect::exact_pcurve_of(&ring_curve, &cap_surface, tol) else {
            ogeom_bail!(Construction, "a cap edge has no closed-form pcurve");
        };
        ogeom_algo::attach_pcurve(
            model,
            &ring0,
            pcurve,
            id,
            ogeom_topo::Location::identity(),
            u_dom,
        )?;
        face
    };

    let faces = [wall, cap];
    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the skinned apex solid did not close");
    }
    make_solid(model, std::slice::from_ref(&sewn.shells[0]))
}

/// A solid skinned over a grid of sections that loops back on itself: the
/// wall is one face closed in both chart directions, no caps at all.
///
/// The `u` seam closes the way every skin's does (pinned row ends), and
/// the `v` loop closes through [`ogeom_geom::fit::fit_surface_grid_closed_v`],
/// C1 across the join. All four boundary traversals are two seam edges used
/// twice, anchored at one shared vertex, exactly as a torus bounds itself.
fn closed_skinned_solid(
    model: &mut Model,
    rows: &[Vec<Point>],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    let shell = closed_skinned_shell(model, rows, tolerance, tol)?;
    make_solid(model, std::slice::from_ref(&shell))
}

/// The closed skin as a shell, for callers assembling solids with voids;
/// a holed profile's ring is one outer shell and one per tunnel.
fn closed_skinned_shell(
    model: &mut Model,
    rows: &[Vec<Point>],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    use ogeom_geom::Surface as _;
    let mut closed_rows: Vec<Vec<Point>> = Vec::with_capacity(rows.len() + 1);
    for row in rows {
        let mut r = row.clone();
        r.push(row[0]);
        closed_rows.push(r);
    }
    closed_rows.push(closed_rows[0].clone());
    let fitted = ogeom_geom::fit::fit_surface_grid_closed_v(&closed_rows, 3, tolerance, tol)?;
    if !fitted.met {
        ogeom_bail!(
            NotDone,
            "the closed skin reached {} against a target of {tolerance}",
            fitted.error
        );
    }
    let surface = fitted.curve;
    let (u_knots, v_knots) = (surface.u_knots().clone(), surface.v_knots().clone());
    let (k, l, net) = {
        let grid = surface.grid();
        let net: Vec<Point> = grid.points().iter().map(|w| (*w).point()).collect();
        (grid.u_count(), grid.v_count(), net)
    };
    let point_at = |i: usize, j: usize| -> Point { net[i * l + j] };
    let (u_dom, v_dom) = surface.domain();

    // Both seams straight off the net: the u-run at v's join, and the v-run
    // at u's.
    let along_u = {
        let control: Vec<Point> = (0..k).map(|i| point_at(i, 0)).collect();
        ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(u_knots, control, tol)?)
    };
    let along_v = {
        let control: Vec<Point> = (0..l).map(|j| point_at(0, j)).collect();
        ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(v_knots, control, tol)?)
    };
    let surface_geo: SurfaceGeometry = surface.into();
    let surface_id = model.geometry_mut().add_surface(surface_geo.clone());

    let u_edge = make_edge(model, along_u, u_dom, tol)?.shape;
    let Some((corner, _)) = ogeom_algo::edge_vertices(model, &u_edge)? else {
        ogeom_bail!(Construction, "the closed skin's seam has no vertex");
    };
    let v_edge = make_edge_between(model, along_v, v_dom, &corner, &corner, tol)?.shape;

    let row_line = |v: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(0.0, v), ogeom_math::Direction2::X),
            u_dom.0 - 1.0,
            u_dom.1 + 1.0,
        )?
        .into())
    };
    let column_line = |u: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(u, 0.0), ogeom_math::Direction2::Y),
            v_dom.0 - 1.0,
            v_dom.1 + 1.0,
        )?
        .into())
    };
    // The u-run is a seam in v (the same curve at both rows), and the
    // v-run a seam in u.
    ogeom_algo::attach_seam(
        model,
        &u_edge,
        row_line(v_dom.0)?,
        row_line(v_dom.1)?,
        surface_id,
        ogeom_topo::Location::identity(),
        u_dom,
    )?;
    ogeom_algo::attach_seam(
        model,
        &v_edge,
        column_line(u_dom.1)?,
        column_line(u_dom.0)?,
        surface_id,
        ogeom_topo::Location::identity(),
        v_dom,
    )?;

    let wire = ogeom_algo::make_wire(
        model,
        &[
            u_edge.clone(),
            v_edge.clone(),
            u_edge.reversed(),
            v_edge.reversed(),
        ],
        tol,
    )?
    .shape;
    let face = ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape;
    // Away from the centre of the section the point lies on; a closed
    // skin's centre as a whole sits in its hole.
    let mid_u = f64::midpoint(u_dom.0, u_dom.1);
    let mid_v = f64::midpoint(v_dom.0, v_dom.1);
    let s_mid = surface_geo.point_at(mid_u, mid_v, tol)?;
    let (du, dv) = surface_geo.d1_at(mid_u, mid_v, tol)?;
    let face = if du.cross(dv).dot(s_mid - section_centre(rows, s_mid)) >= 0.0 {
        face
    } else {
        face.reversed()
    };

    let sewn = sew(model, std::slice::from_ref(&face), tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the closed skin did not close");
    }
    Ok(sewn.shells[0].clone())
}

/// The loft to a point: a wire section closing onto a single apex vertex.
///
/// A circle takes the cone the revolved primitives already build, apex on
/// its axis or refused; a polygon takes exact planar triangle walls, sound
/// for *any* apex: a skew pyramid's walls are still triangles.
fn loft_to_point(
    model: &mut Model,
    section: &Shape,
    apex: &Shape,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if !ogeom_algo::is_wire_closed(model, section, tol)? {
        ogeom_bail!(Construction, "a loft section must be closed");
    }
    let apex_point = {
        let Some(data) = model.node(apex).and_then(|n| n.data().as_vertex()) else {
            ogeom_bail!(Construction, "the apex vertex holds no data");
        };
        data.point
    };

    // The circular case: a cone, apex on the axis.
    let edges = explore(model, section, Filter::OfType(ShapeType::Edge))?;
    if edges.len() == 1
        && let Some(data) = model.node(&edges[0]).and_then(|n| n.data().as_edge())
        && let Some(EdgeRepr::Curve3d { curve, .. }) = data.curve3d()
        && let Some(Curve::Circle(c)) = model.geometry().curve(*curve)
    {
        let circle = c.circle();
        let axis = circle.frame().z().vector();
        let rise = apex_point - circle.centre();
        let height = rise.dot(axis);
        if rise.cross(axis).magnitude() > tol.confusion() * 10.0 {
            ogeom_bail!(
                Construction,
                "a circle lofts to a point on its own axis; the oblique cone \
                 needs the skinned machinery; see docs/PARITY.md, offset.loft"
            );
        }
        if height.abs() <= tol.confusion() {
            ogeom_bail!(Construction, "the apex sits in the section's own plane");
        }
        let base = if height > 0.0 {
            circle.frame()
        } else {
            Frame::new(
                circle.centre(),
                -circle.frame().z(),
                circle.frame().x(),
                tol,
            )?
        };
        let mut built =
            ogeom_algo::make_cone(model, base, circle.radius(), 0.0, height.abs(), tol)?;
        built.history.generate(section, built.shape.clone());
        built.history.generate(apex, built.shape.clone());
        return Ok(built);
    }

    // The polygonal case: exact triangle walls to a shared apex.
    let mut corners: Vec<Point> = Vec::new();
    for edge in model.ordered_children_of(section)? {
        let Some(data) = model.node(&edge).and_then(|n| n.data().as_edge()) else {
            ogeom_bail!(Construction, "a section edge holds no data");
        };
        let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
            ogeom_bail!(Construction, "a section edge has no curve");
        };
        let Some(Curve::Line(line)) = model.geometry().curve(*curve).cloned() else {
            ogeom_bail!(
                Construction,
                "a mixed or curved section lofts to a point through the \
                 skinned machinery; see docs/PARITY.md, offset.loft"
            );
        };
        let t = if edge.orientation() == ogeom_topo::Orientation::Reversed {
            range.1
        } else {
            range.0
        };
        corners.push(ogeom_geom::Curve::Line(line).point_at(t, tol)?);
    }
    if corners.len() < 3 {
        ogeom_bail!(Construction, "a pyramid needs at least three base corners");
    }
    let apex_vertex = ogeom_algo::make_vertex(model, apex_point).shape;
    let base_vertices: Vec<Shape> = corners
        .iter()
        .map(|p| ogeom_algo::make_vertex(model, *p).shape)
        .collect();
    let segment =
        |model: &mut Model, from: (&Shape, Point), to: (&Shape, Point)| -> OgeomResult<Shape> {
            let line = ogeom_geom::LineCurve::segment(from.1, to.1, tol)?;
            let curve: Curve = line.into();
            let domain = curve.domain();
            Ok(make_edge_between(model, curve, domain, from.0, to.0, tol)?.shape)
        };
    let count = corners.len();
    let mut base_edges = Vec::with_capacity(count);
    let mut rails = Vec::with_capacity(count);
    for i in 0..count {
        let next = (i + 1) % count;
        base_edges.push(segment(
            model,
            (&base_vertices[i], corners[i]),
            (&base_vertices[next], corners[next]),
        )?);
        rails.push(segment(
            model,
            (&base_vertices[i], corners[i]),
            (&apex_vertex, apex_point),
        )?);
    }
    let centroid = {
        let mut c = Vector::new(0.0, 0.0, 0.0);
        for p in &corners {
            c += p.to_vector();
        }
        #[allow(clippy::cast_precision_loss)]
        Point::from_vector(c / count as f64 / 4.0 * 3.0 + apex_point.to_vector() / 4.0)
    };
    let planar = |model: &mut Model, pts: [Point; 3], walk: Vec<Shape>| -> OgeomResult<Shape> {
        let n = (pts[1] - pts[0]).cross(pts[2] - pts[0]);
        let m = n.magnitude();
        if m <= tol.confusion() {
            ogeom_bail!(Construction, "a wall of the pyramid is degenerate");
        }
        let mut outward = n / m;
        if outward.dot(pts[0] - centroid) < 0.0 {
            outward = -outward;
        }
        let plane = ogeom_math::Plane::through(pts[0], Direction::new(outward, tol)?);
        let mut reach = 1.0_f64;
        for p in pts {
            reach = reach.max(p.distance(pts[0]) * 2.0);
        }
        let surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        let id = model.geometry_mut().add_surface(surface.clone());
        let signed = {
            let (du, dv) = {
                use ogeom_geom::Surface as _;
                surface.d1_at(0.0, 0.0, tol)?
            };
            du.cross(dv).dot(outward) >= 0.0
        };
        let mut wired = Vec::with_capacity(walk.len());
        for used in &walk {
            let (curve, range) = spine_curve_of(model, used)?;
            let Some(pcurve) = ogeom_intersect::exact_pcurve_of(&curve, &surface, tol) else {
                ogeom_bail!(Construction, "a wall edge has no closed-form pcurve");
            };
            ogeom_algo::attach_pcurve(
                model,
                used,
                pcurve,
                id,
                ogeom_topo::Location::identity(),
                range,
            )?;
            wired.push(used.clone());
        }
        let wire = ogeom_algo::make_wire(model, &wired, tol)?.shape;
        let face = ogeom_algo::make_face_on(model, id, std::slice::from_ref(&wire), tol)?.shape;
        Ok(if signed { face } else { face.reversed() })
    };
    let mut faces = Vec::with_capacity(count + 1);
    for i in 0..count {
        let next = (i + 1) % count;
        faces.push(planar(
            model,
            [corners[i], corners[next], apex_point],
            vec![
                base_edges[i].clone(),
                rails[next].clone(),
                rails[i].reversed(),
            ],
        )?);
    }
    // The base cap: all corners, wound against the walls.
    let base_walk: Vec<Shape> = (0..count).rev().map(|i| base_edges[i].reversed()).collect();
    faces.push({
        let n = (corners[1] - corners[0]).cross(corners[2] - corners[0]);
        let mut outward = n / n.magnitude();
        if outward.dot(corners[0] - centroid) < 0.0 {
            outward = -outward;
        }
        let plane = ogeom_math::Plane::through(corners[0], Direction::new(outward, tol)?);
        let mut reach = 1.0_f64;
        for p in &corners {
            reach = reach.max(p.distance(corners[0]) * 2.0);
        }
        let surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        let id = model.geometry_mut().add_surface(surface.clone());
        for used in &base_walk {
            let (curve, range) = spine_curve_of(model, used)?;
            let Some(pcurve) = ogeom_intersect::exact_pcurve_of(&curve, &surface, tol) else {
                ogeom_bail!(Construction, "a base edge has no closed-form pcurve");
            };
            ogeom_algo::attach_pcurve(
                model,
                used,
                pcurve,
                id,
                ogeom_topo::Location::identity(),
                range,
            )?;
        }
        let wire = ogeom_algo::make_wire(model, &base_walk, tol)?.shape;
        let face = ogeom_algo::make_face_on(model, id, std::slice::from_ref(&wire), tol)?.shape;
        let signed = {
            use ogeom_geom::Surface as _;
            let (du, dv) = surface.d1_at(0.0, 0.0, tol)?;
            du.cross(dv).dot(outward) >= 0.0
        };
        if signed { face } else { face.reversed() }
    });

    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the pyramid did not close");
    }
    let mut built = make_solid(model, std::slice::from_ref(&sewn.shells[0]))?;
    built.history.generate(section, built.shape.clone());
    built.history.generate(apex, built.shape.clone());
    Ok(built)
}

/// Loft through sections with the start of each row named by the caller.
///
/// [`make_loft_skinned`] leaves alignment to each section's own traversal
/// start; this sibling takes one hint per section (a point near where its
/// row should begin) and rotates each sampling there, which is how a
/// caller untwists a loft whose wires happen to start in different places.
///
/// # Errors
///
/// As [`make_loft_skinned`], and additionally if the hints do not pair up
/// with the sections.
pub fn make_loft_skinned_aligned(
    model: &mut Model,
    sections: &[Shape],
    hints: &[Point],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if hints.len() != sections.len() {
        ogeom_bail!(
            Construction,
            "{} hints against {} sections; each section names its own start",
            hints.len(),
            sections.len()
        );
    }
    if sections.len() < 2 {
        ogeom_bail!(Construction, "a loft needs at least two sections");
    }
    const AROUND: usize = 48;
    let mut rows: Vec<Vec<Point>> = Vec::with_capacity(sections.len());
    let mut planes: Vec<Plane> = Vec::with_capacity(sections.len());
    for (wire, hint) in sections.iter().zip(hints) {
        if model.kind_of(wire)? != ShapeType::Wire {
            ogeom_bail!(Construction, "a loft section is a closed wire");
        }
        if !ogeom_algo::is_wire_closed(model, wire, tol)? {
            ogeom_bail!(Construction, "a loft section must be closed");
        }
        let Some(plane) = ogeom_algo::find_plane(model, wire, tol)? else {
            ogeom_bail!(Construction, "a loft section must be planar");
        };
        planes.push(plane);
        rows.push(sample_wire_from(model, wire, AROUND, Some(*hint), tol)?);
    }
    let outward0 = {
        let towards = rows[1][0] - rows[0][0];
        let n = planes[0].normal().vector();
        if n.dot(towards) > 0.0 { -n } else { n }
    };
    let outward1 = {
        let towards = rows[rows.len() - 2][0] - rows[rows.len() - 1][0];
        let n = planes[planes.len() - 1].normal().vector();
        if n.dot(towards) > 0.0 { -n } else { n }
    };
    let mut built = skinned_solid(
        model,
        &rows,
        (EndCap::Plane(outward0), EndCap::Plane(outward1)),
        tolerance,
        tol,
    )?;
    for section in sections {
        built.history.generate(section, built.shape.clone());
    }
    Ok(built)
}

/// Loft a ring through closed planar sections that loop back to the first.
///
/// [`make_loft_skinned`]'s closed sibling: the sections are sampled the same
/// way, the skin runs through all of them and back to the start, C1 across
/// the loop, and there are no caps: the result bounds itself the way a
/// torus does. The sections are *not* repeated: the loop-back is the
/// construction's own.
///
/// The closed join costs freedom: a sparse loop fits only loosely, and the
/// refusal quotes the deviation it honestly reached. A loop that wants a
/// tight tolerance wants sections dense enough to bend around: in
/// practice, a dozen and up.
///
/// # Errors
///
/// As [`make_loft_skinned`], needing at least three sections;
/// [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) if the closed
/// skin cannot reach the tolerance.
pub fn make_loft_skinned_closed(
    model: &mut Model,
    sections: &[Shape],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if sections.len() < 3 {
        ogeom_bail!(Construction, "a closed loft needs at least three sections");
    }
    const AROUND: usize = 48;
    let mut rows: Vec<Vec<Point>> = Vec::with_capacity(sections.len());
    for wire in sections {
        if model.kind_of(wire)? != ShapeType::Wire {
            ogeom_bail!(Construction, "a loft section is a closed wire");
        }
        if !ogeom_algo::is_wire_closed(model, wire, tol)? {
            ogeom_bail!(Construction, "a loft section must be closed");
        }
        if ogeom_algo::find_plane(model, wire, tol)?.is_none() {
            ogeom_bail!(Construction, "a loft section must be planar");
        }
        rows.push(sample_wire(model, wire, AROUND, tol)?);
    }
    let mut built = closed_skinned_solid(model, &rows, tolerance, tol)?;
    for section in sections {
        built.history.generate(section, built.shape.clone());
    }
    Ok(built)
}

/// A skinned strip: one open patch of a sweep, with its border edges.
///
/// The wall of a *faceted* profile cannot be one closed skin (a fit cannot
/// speak a corner), so each profile edge sweeps its own strip, cornered at
/// the caller's shared vertices, and the strips weld along their rails by
/// the tolerance the fit honestly carries.
struct SkinnedStrip {
    face: Shape,
    /// The border along the first station, from `corners.0` to `corners.1`.
    bottom: Shape,
    /// The border along the last station, from `corners.2` to `corners.3`.
    top: Shape,
    /// The rail along the profile edge's start, from `corners.0` to `corners.2`.
    rail0: Shape,
    /// The rail along the profile edge's end, from `corners.1` to `corners.3`.
    rail1: Shape,
}

/// Skin an open grid of samples (stations by profile-edge samples) into
/// one strip. `corners` are the caller's vertices at (first station, edge
/// start), (first, end), (last, start), (last, end), shared with the
/// neighbouring strips so the wires chain. `shared` are borders a
/// neighbouring strip already built (bottom, top, start rail, end rail),
/// adopted as they are (see `adopt_border`).
#[allow(clippy::too_many_arguments, reason = "one strip, spelled out")]
fn skinned_strip(
    model: &mut Model,
    rows: &[Vec<Point>],
    corners: (&Shape, &Shape, &Shape, &Shape),
    shared: [Option<&Shape>; 4],
    outward_hint: Point,
    hole: bool,
    by_spacing: bool,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<SkinnedStrip> {
    use ogeom_geom::Surface as _;
    // A strip whose every row lies in one plane (a straight profile edge
    // down a straight run, a flat face of the profile along a planar
    // spine) is that plane, exactly: a coplanar neighbour then melts with
    // it on the one surface two fits of it would never agree on.
    if let Some(plane) = plane_of_rows(rows, tol) {
        return planar_strip(
            model,
            rows,
            plane,
            corners,
            shared,
            outward_hint,
            hole,
            tolerance,
            tol,
        );
    }
    // Sections a caller placed follow their own spacing across the skin;
    // stations a sweep placed keep the fit's centripetal assignment.
    let fitted = if by_spacing {
        ogeom_geom::fit::fit_surface_grid_sections(rows, 3, tolerance, tol)?
    } else {
        ogeom_geom::fit::fit_surface_grid(rows, 3, tolerance, tol)?
    };
    if !fitted.met {
        ogeom_bail!(
            NotDone,
            "the strip reached {} against a target of {tolerance}",
            fitted.error
        );
    }
    let error = fitted.error.max(tol.confusion());
    let surface = fitted.curve;
    let (u_knots, v_knots) = (surface.u_knots().clone(), surface.v_knots().clone());
    let (k, l, net) = {
        let grid = surface.grid();
        let net: Vec<Point> = grid.points().iter().map(|w| (*w).point()).collect();
        (grid.u_count(), grid.v_count(), net)
    };
    let point_at = |i: usize, j: usize| -> Point { net[i * l + j] };
    let (u_dom, v_dom) = surface.domain();

    let u_curve = |j: usize| -> OgeomResult<ogeom_geom::Curve> {
        let control: Vec<Point> = (0..k).map(|i| point_at(i, j)).collect();
        Ok(ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            u_knots.clone(),
            control,
            tol,
        )?))
    };
    let v_curve = |i: usize| -> OgeomResult<ogeom_geom::Curve> {
        let control: Vec<Point> = (0..l).map(|j| point_at(i, j)).collect();
        Ok(ogeom_geom::Curve::BSpline(ogeom_geom::BSplineCurve::new(
            v_knots.clone(),
            control,
            tol,
        )?))
    };
    let surface_geo: SurfaceGeometry = surface.into();
    let surface_id = model.geometry_mut().add_surface(surface_geo.clone());

    let (c00, c10, c01, c11) = corners;
    // A border the neighbouring strip already built is adopted, not
    // refitted (see `adopt_border`).
    let border = |model: &mut Model,
                  given: Option<&Shape>,
                  curve: ogeom_geom::Curve,
                  range: (f64, f64),
                  from: &Shape,
                  to: &Shape|
     -> OgeomResult<Shape> {
        match given {
            Some(edge) => {
                adopt_border(model, edge, &surface_geo, error + tol.confusion(), tol)?;
                Ok(edge.clone())
            }
            None => Ok(make_edge_between(model, curve, range, from, to, tol)?.shape),
        }
    };
    let bottom = border(model, shared[0], u_curve(0)?, u_dom, c00, c10)?;
    let top = border(model, shared[1], u_curve(l - 1)?, u_dom, c01, c11)?;
    let rail0 = border(model, shared[2], v_curve(0)?, v_dom, c00, c01)?;
    let rail1 = border(model, shared[3], v_curve(k - 1)?, v_dom, c10, c11)?;

    let row_line = |v: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(0.0, v), ogeom_math::Direction2::X),
            u_dom.0 - 1.0,
            u_dom.1 + 1.0,
        )?
        .into())
    };
    let column_line = |u: f64| -> OgeomResult<ogeom_geom::PlanarCurve> {
        Ok(Line2d::over(
            ogeom_math::Axis2::new(Point2::new(u, 0.0), ogeom_math::Direction2::Y),
            v_dom.0 - 1.0,
            v_dom.1 + 1.0,
        )?
        .into())
    };
    // A border this strip fitted runs along its own row or column, and its
    // image is that straight line. One adopted from the neighbour was fitted
    // at the neighbour's pace along the sweep, which is not this strip's:
    // its image here is read off this surface point by point.
    for (edge, given, straight, span) in [
        (&bottom, shared[0].is_some(), row_line(v_dom.0)?, u_dom),
        (&top, shared[1].is_some(), row_line(v_dom.1)?, u_dom),
        (&rail0, shared[2].is_some(), column_line(u_dom.0)?, v_dom),
        (&rail1, shared[3].is_some(), column_line(u_dom.1)?, v_dom),
    ] {
        let (image, range) = if given {
            adopted_image(model, edge, &surface_geo, tol)?
        } else {
            (straight, span)
        };
        ogeom_algo::attach_pcurve(
            model,
            edge,
            image,
            surface_id,
            ogeom_topo::Location::identity(),
            range,
        )?;
    }
    // The rails carry the fit's honest budget: the neighbouring strip fitted
    // the same transported corners independently, and the weld between them
    // is only as tight as both fits.
    for edge in [&bottom, &top, &rail0, &rail1] {
        model.widen(edge, ogeom_core::Tolerance::new(error)?)?;
    }

    let wire = ogeom_algo::make_wire(
        model,
        &[
            bottom.clone(),
            rail1.clone(),
            top.reversed(),
            rail0.reversed(),
        ],
        tol,
    )?
    .shape;
    let face = ogeom_algo::make_face_on(model, surface_id, std::slice::from_ref(&wire), tol)?.shape;
    let mid_u = f64::midpoint(u_dom.0, u_dom.1);
    let mid_v = f64::midpoint(v_dom.0, v_dom.1);
    let s_mid = surface_geo.point_at(mid_u, mid_v, tol)?;
    let (du, dv) = surface_geo.d1_at(mid_u, mid_v, tol)?;
    let natural_out = du.cross(dv).dot(s_mid - outward_hint) >= 0.0;
    let face = if natural_out == !hole {
        face
    } else {
        face.reversed()
    };
    Ok(SkinnedStrip {
        face,
        bottom,
        top,
        rail0,
        rail1,
    })
}

/// The plane every point of the rows lies in, if there is one.
fn plane_of_rows(rows: &[Vec<Point>], tol: Tolerances) -> Option<Plane> {
    let first = rows.first()?;
    let last = rows.last()?;
    let origin = *first.first()?;
    let across = *first.last()? - origin;
    let along = *last.first()? - origin;
    let normal = across.cross(along);
    if normal.magnitude() <= tol.confusion() * across.magnitude().max(along.magnitude()) {
        return None;
    }
    let normal = Direction::new(normal, tol).ok()?;
    let plane = Plane::through(origin, normal);
    rows.iter()
        .flatten()
        .all(|p| plane.distance_to(*p) <= tol.confusion())
        .then_some(plane)
}

/// A strip on its own exact plane: the borders fitted through the rows and
/// the end columns, the face on the plane.
#[allow(clippy::too_many_arguments, reason = "one construction, all its data")]
fn planar_strip(
    model: &mut Model,
    rows: &[Vec<Point>],
    plane: Plane,
    corners: (&Shape, &Shape, &Shape, &Shape),
    shared: [Option<&Shape>; 4],
    outward_hint: Point,
    hole: bool,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<SkinnedStrip> {
    // Splines, as every swept border is: the caps and the neighbouring
    // strips read them so, and a spline through collinear points is the
    // straight segment itself.
    let through = |points: &[Point]| -> OgeomResult<ogeom_geom::Curve> {
        let fitted = ogeom_geom::fit::fit_points(points, 3, tolerance * 0.5, tol)?;
        if !fitted.met {
            ogeom_bail!(
                NotDone,
                "a planar strip's border reached {} against a target of {tolerance}",
                fitted.error
            );
        }
        Ok(ogeom_geom::Curve::BSpline(fitted.curve))
    };
    let column = |i: usize| -> Vec<Point> { rows.iter().map(|row| row[i]).collect() };
    let last = rows[0].len() - 1;
    let (c00, c10, c01, c11) = corners;
    let border = |model: &mut Model,
                  given: Option<&Shape>,
                  points: Vec<Point>,
                  from: &Shape,
                  to: &Shape|
     -> OgeomResult<Shape> {
        if let Some(edge) = given {
            return Ok(edge.clone());
        }
        let curve = through(&points)?;
        let domain = curve.domain();
        Ok(make_edge_between(model, curve, domain, from, to, tol)?.shape)
    };
    let bottom = border(model, shared[0], rows[0].clone(), c00, c10)?;
    let top = border(model, shared[1], rows[rows.len() - 1].clone(), c01, c11)?;
    let rail0 = border(model, shared[2], column(0), c00, c01)?;
    let rail1 = border(model, shared[3], column(last), c10, c11)?;

    // The loop runs across, up, back and down: counter-clockwise about the
    // normal from the first row's run to the first column's.
    let across = rows[0][last] - rows[0][0];
    let up = rows[rows.len() - 1][0] - rows[0][0];
    let normal = Direction::new(across.cross(up), tol)?;
    let wound = Plane::through(plane.origin(), normal);
    let reach = rows
        .iter()
        .flatten()
        .map(|p| p.distance(plane.origin()))
        .fold(1.0_f64, f64::max)
        * 2.0;
    let surface: SurfaceGeometry =
        PlaneSurface::over(wound, (-reach, reach), (-reach, reach))?.into();
    let face = ogeom_algo::make_face_with_pcurves(
        model,
        surface.clone(),
        &[vec![
            bottom.clone(),
            rail1.clone(),
            top.reversed(),
            rail0.reversed(),
        ]],
        tol,
    )?
    .shape;
    let mid = rows[rows.len() / 2][last / 2];
    let natural_out = normal.vector().dot(mid - outward_hint) >= 0.0;
    let face = if natural_out == !hole {
        face
    } else {
        face.reversed()
    };
    Ok(SkinnedStrip {
        face,
        bottom,
        top,
        rail0,
        rail1,
    })
}

/// A loft through circles standing coaxial on parallel planes: the solid of
/// revolution of the meridian through their radii, a spline through them
/// in the half-plane of the first circle's start. `None` where the
/// sections are anything else.
fn coaxial_circles_loft(
    model: &mut Model,
    sections: &[Shape],
    tol: Tolerances,
) -> OgeomResult<Option<Built>> {
    let mut circles = Vec::with_capacity(sections.len());
    for wire in sections {
        if model.kind_of(wire)? != ShapeType::Wire {
            return Ok(None);
        }
        let edges = model.ordered_children_of(wire)?;
        let [edge] = edges.as_slice() else {
            return Ok(None);
        };
        let (curve, _) = spine_curve_of(model, edge)?;
        let ogeom_geom::Curve::Circle(c) = curve else {
            return Ok(None);
        };
        let placed = c
            .circle()
            .transformed(&edge.transform(model.datums())?, tol)?;
        circles.push(placed);
    }
    let first = circles[0].frame();
    let (c0, z0) = (first.origin(), first.z().vector());
    let last = circles[circles.len() - 1].centre();
    let rise = last - c0;
    if rise.magnitude() <= tol.confusion() {
        return Ok(None);
    }
    let z = rise / rise.magnitude();
    if z.cross(z0).magnitude() > tol.angular() {
        return Ok(None);
    }
    let mut heights = Vec::with_capacity(circles.len());
    for c in &circles {
        let off = c.centre() - c0;
        if off.cross(z).magnitude() > tol.confusion() * 10.0
            || c.frame().z().vector().cross(z).magnitude() > tol.angular()
        {
            return Ok(None);
        }
        heights.push(off.dot(z));
    }
    if heights.windows(2).any(|w| w[1] <= w[0] + tol.confusion()) {
        return Ok(None);
    }
    let x = first.x().vector();
    let meridian: Vec<Point> = circles
        .iter()
        .zip(&heights)
        .map(|(c, h)| c0 + z * *h + x * c.radius())
        .collect();
    let degree = (meridian.len() - 1).min(3);
    let fitted = ogeom_geom::fit::fit_points(&meridian, degree, tol.confusion() * 1e-3, tol)?;
    let spline: ogeom_geom::Curve = fitted.curve.into();
    let domain = spline.domain();
    let top = c0 + z * heights[heights.len() - 1];
    let vertex = |model: &mut Model, p: Point| ogeom_algo::make_vertex(model, p).shape;
    let (v_axis0, v_axis1) = (vertex(model, c0), vertex(model, top));
    let (v_rim0, v_rim1) = (
        vertex(model, meridian[0]),
        vertex(model, meridian[meridian.len() - 1]),
    );
    let segment =
        |model: &mut Model, a: (&Shape, Point), b: (&Shape, Point)| -> OgeomResult<Shape> {
            let line: ogeom_geom::Curve = LineCurve::segment(a.1, b.1, tol)?.into();
            let range = line.domain();
            Ok(make_edge_between(model, line, range, a.0, b.0, tol)?.shape)
        };
    let bottom = segment(model, (&v_axis0, c0), (&v_rim0, meridian[0]))?;
    let side = make_edge_between(model, spline, domain, &v_rim0, &v_rim1, tol)?.shape;
    let top_edge = segment(
        model,
        (&v_rim1, meridian[meridian.len() - 1]),
        (&v_axis1, top),
    )?;
    let axis_edge = segment(model, (&v_axis1, top), (&v_axis0, c0))?;
    let wire = ogeom_algo::make_wire(model, &[bottom, side, top_edge, axis_edge], tol)?.shape;
    // Framed from a point inside the profile: a plane's own origin is
    // where its face is read when it carries no trims.
    let inside = c0
        + z * (heights[heights.len() - 1] * 0.5)
        + x * (circles
            .iter()
            .map(|c| c.radius())
            .fold(f64::INFINITY, f64::min)
            * 0.5);
    let plane = Plane::new(Frame::new(
        inside,
        Direction::new(z.cross(x), tol)?,
        Direction::new(x, tol)?,
        tol,
    )?);
    let face = ogeom_algo::make_face(model, PlaneSurface::new(plane).into(), &[wire], tol)?.shape;
    let axis = ogeom_math::Axis {
        location: c0,
        direction: Direction::new(z, tol)?,
    };
    let built = ogeom_algo::make_revolution(model, &face, axis, core::f64::consts::TAU, tol)?;
    Ok(Some(built))
}

/// A loft through sections of one edge count, every vertex a corner: one
/// strip per edge through all the sections, meeting its neighbours along
/// seams through the matched corners, each strip a plane wherever its rows
/// share one. `None` where the sections do not pair edge for edge.
fn cornered_loft(
    model: &mut Model,
    sections: &[Shape],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Option<Built>> {
    let mut rings: Vec<Vec<Shape>> = Vec::with_capacity(sections.len());
    for wire in sections {
        if model.kind_of(wire)? != ShapeType::Wire || !ogeom_algo::is_wire_closed(model, wire, tol)?
        {
            return Ok(None);
        }
        rings.push(model.ordered_children_of(wire)?);
    }
    let count = rings[0].len();
    if count < 2 || rings.iter().any(|r| r.len() != count) {
        return Ok(None);
    }
    let (Some(plane0), Some(plane1)) = (
        ogeom_algo::find_plane(model, &sections[0], tol)?,
        ogeom_algo::find_plane(model, &sections[sections.len() - 1], tol)?,
    ) else {
        return Ok(None);
    };
    const ALONG: usize = 16;
    // Every edge of every section sampled in its ring's sense.
    let mut samples: Vec<Vec<Vec<Point>>> = Vec::with_capacity(rings.len());
    for ring in &rings {
        let mut per_edge = Vec::with_capacity(count);
        for edge in ring {
            let (curve, range) = spine_curve_of(model, edge)?;
            let curve = curve.transformed(&edge.transform(model.datums())?, tol)?;
            let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
            let mut row = Vec::with_capacity(ALONG + 1);
            for k in 0..=ALONG {
                #[allow(clippy::cast_precision_loss)]
                let f = k as f64 / ALONG as f64;
                let t = if reversed {
                    range.1 - (range.1 - range.0) * f
                } else {
                    range.0 + (range.1 - range.0) * f
                };
                row.push(curve.point_at(t, tol)?);
            }
            per_edge.push(row);
        }
        samples.push(per_edge);
    }
    let corners = |model: &mut Model, s: usize| -> Vec<Shape> {
        (0..count)
            .map(|e| ogeom_algo::make_vertex(model, samples[s][e][0]).shape)
            .collect()
    };
    let (from, to) = (corners(model, 0), corners(model, sections.len() - 1));
    let middle = &samples[sections.len() / 2];
    let hint = {
        let all: Vec<Point> = middle.iter().flatten().copied().collect();
        #[allow(clippy::cast_precision_loss)]
        let n = all.len() as f64;
        Point::from_vector(
            all.iter()
                .fold(Vector::new(0.0, 0.0, 0.0), |acc, p| acc + p.to_vector())
                / n,
        )
    };
    let mut faces = Vec::with_capacity(count + 2);
    let (mut bottoms, mut tops) = (Vec::with_capacity(count), Vec::with_capacity(count));
    let mut first_rail: Option<Shape> = None;
    let mut prev_rail: Option<Shape> = None;
    for e in 0..count {
        let rows: Vec<Vec<Point>> = samples.iter().map(|s| s[e].clone()).collect();
        let next = (e + 1) % count;
        let last_rail = if e + 1 == count {
            first_rail.clone()
        } else {
            None
        };
        let strip = skinned_strip(
            model,
            &rows,
            (&from[e], &from[next], &to[e], &to[next]),
            [None, None, prev_rail.as_ref(), last_rail.as_ref()],
            hint,
            false,
            true,
            tolerance,
            tol,
        )?;
        if e == 0 {
            first_rail = Some(strip.rail0.clone());
        }
        prev_rail = Some(strip.rail1.clone());
        faces.push(strip.face.clone());
        bottoms.push(strip.bottom);
        tops.push(strip.top);
    }
    let towards = hint - samples[0][0][0];
    let n0 = plane0.normal().vector();
    let n0 = if n0.dot(towards) > 0.0 { -n0 } else { n0 };
    let away = hint - samples[sections.len() - 1][0][0];
    let n1 = plane1.normal().vector();
    let n1 = if n1.dot(away) > 0.0 { -n1 } else { n1 };
    faces.push(plane_cap(model, samples[0][0][0], n0, &[bottoms], tol)?);
    faces.push(plane_cap(
        model,
        samples[sections.len() - 1][0][0],
        n1,
        &[tops],
        tol,
    )?);
    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the cornered loft did not close");
    }
    Ok(Some(make_solid(model, &sewn.shells)?))
}

/// A planar cap through `at`, facing `outward`, bounded by loops of spline
/// or line edges lying in it; each edge's trim is its exact projection into
/// the plane's chart.
fn plane_cap(
    model: &mut Model,
    at: Point,
    outward: Vector,
    loops: &[Vec<Shape>],
    tol: Tolerances,
) -> OgeomResult<Shape> {
    let cap_plane = Plane::through(at, Direction::new(outward, tol)?);
    let mut reach = 1.0_f64;
    for edges in loops {
        for edge in edges {
            let (curve, range) = spine_curve_of(model, edge)?;
            for k in 0..8 {
                let p = curve.point_at(range.0 + (range.1 - range.0) * f64::from(k) / 8.0, tol)?;
                reach = reach.max(p.distance(at) * 2.0);
            }
        }
    }
    let surface: SurfaceGeometry =
        PlaneSurface::over(cap_plane, (-reach, reach), (-reach, reach))?.into();
    let mut wires = Vec::with_capacity(loops.len());
    for edges in loops {
        wires.push(ogeom_algo::make_wire(model, edges, tol)?.shape);
    }
    let face = ogeom_algo::make_face(model, surface, &wires, tol)?.shape;
    let cap_id = {
        let Some(ogeom_topo::NodeData::Face(data)) = model.node(&face).map(|n| n.data()) else {
            ogeom_bail!(Construction, "the cap holds no face data");
        };
        data.surface
    };
    let frame = cap_plane.frame();
    let flat = |p: Point| {
        let local = frame.to_local(p);
        Point2::new(local.x, local.y)
    };
    for edges in loops {
        for edge in edges {
            let (curve, range) = spine_curve_of(model, edge)?;
            let pcurve: ogeom_geom::PlanarCurve = match &curve {
                ogeom_geom::Curve::BSpline(bs) => {
                    let control2: Vec<Point2> = bs
                        .control_points()
                        .iter()
                        .map(|w| flat(w.point()))
                        .collect();
                    ogeom_geom::BSpline2d::new(bs.knots().clone(), control2, tol)?.into()
                }
                ogeom_geom::Curve::Line(line) => {
                    let axis = line.axis();
                    let origin = flat(axis.location);
                    let ahead = flat(axis.location + axis.direction.vector());
                    ogeom_geom::Line2d::over(
                        ogeom_math::Axis2::through(origin, ahead, tol)?,
                        range.0,
                        range.1,
                    )?
                    .into()
                }
                _ => ogeom_bail!(Construction, "a cap edge is neither a spline nor a line"),
            };
            ogeom_algo::attach_pcurve(
                model,
                edge,
                pcurve,
                cap_id,
                ogeom_topo::Location::identity(),
                range,
            )?;
        }
    }
    Ok(face)
}

/// Loft a solid through many closed planar sections, skinned smoothly.
///
/// The sections are sampled at matched arc-length fractions from their own
/// traversal starts (aligning those starts is the caller's authorship),
/// and the skin holds every section to `tolerance`. The caps are the first
/// and last sections' own planes.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if fewer than
/// two sections, a section is open or not planar;
/// [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) if
/// the skin cannot reach the tolerance.
pub fn make_loft_skinned(
    model: &mut Model,
    sections: &[Shape],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if sections.len() < 2 {
        ogeom_bail!(Construction, "a loft needs at least two sections");
    }
    let to_point = model.kind_of(&sections[sections.len() - 1])? == ShapeType::Vertex;
    // The smooth skin through two sections is the ruled one, and that is
    // built exactly: a drum or cone between circles, planes between
    // polygons.
    if sections.len() == 2
        && !to_point
        && let Ok(built) = make_loft(model, &sections[0], &sections[1], tol)
    {
        return Ok(built);
    }
    if sections.len() > 2 && !to_point {
        if let Some(built) = coaxial_circles_loft(model, sections, tol)? {
            return Ok(built);
        }
        if let Some(built) = cornered_loft(model, sections, tolerance, tol)? {
            return Ok(built);
        }
    }
    const AROUND: usize = 48;
    // A trailing vertex is the apex form: the skin narrows to a point and
    // the solid closes there without a cap.
    let apex = match model.kind_of(&sections[sections.len() - 1])? {
        ShapeType::Vertex => {
            if sections.len() < 2 {
                ogeom_bail!(
                    Construction,
                    "a loft to a point needs a section to start from"
                );
            }
            let Some(data) = model
                .node(&sections[sections.len() - 1])
                .and_then(|n| n.data().as_vertex())
            else {
                ogeom_bail!(Construction, "the apex vertex holds no point");
            };
            Some(data.point)
        }
        _ => None,
    };
    let wires = &sections[..sections.len() - usize::from(apex.is_some())];
    let mut rows: Vec<Vec<Point>> = Vec::with_capacity(sections.len());
    let mut cap_planes: Vec<Option<Plane>> = Vec::with_capacity(wires.len());
    for wire in wires {
        if model.kind_of(wire)? != ShapeType::Wire {
            ogeom_bail!(Construction, "a loft section is a wire");
        }
        if !ogeom_algo::is_wire_closed(model, wire, tol)? {
            ogeom_bail!(Construction, "a loft section must be closed");
        }
        // Planarity is a *cap's* requirement, not the fit's: only the
        // sections a cap will stand on must hold a plane. A wavy middle
        // section skins fine.
        cap_planes.push(ogeom_algo::find_plane(model, wire, tol)?);
        rows.push(sample_wire(model, wire, AROUND, tol)?);
    }
    // A planar end is capped by its plane; one that is not (a wavy rim),
    // by a patch skinned from the ring to a point inside it.
    let outward_at = |rows: &[Vec<Point>], planes: &[Option<Plane>], end: bool| -> EndCap {
        let (i, j) = if end {
            (rows.len() - 1, rows.len() - 2)
        } else {
            (0, 1)
        };
        let Some(plane) = &planes[i] else {
            return EndCap::Skinned;
        };
        let towards = rows[j][0] - rows[i][0];
        let n = plane.normal().vector();
        EndCap::Plane(if n.dot(towards) > 0.0 { -n } else { n })
    };
    let mut built = if let Some(apex) = apex {
        if rows.len() < 2 {
            // One ring to a point is exact machinery's job when it can be;
            // the skin still needs two rows to shape the wall, so a middle
            // row is interpolated halfway toward the apex.
            let half: Vec<Point> = rows[0]
                .iter()
                .map(|p| Point::from_vector((p.to_vector() + apex.to_vector()) * 0.5))
                .collect();
            rows.push(half);
        }
        let outward0 = match outward_at(&rows, &cap_planes, false) {
            EndCap::Plane(n) => n,
            EndCap::Skinned => {
                ogeom_bail!(
                    Construction,
                    "a loft to a point starts from a planar section; a cap stands on it"
                );
            }
        };
        rows.push(vec![apex; AROUND]);
        skinned_solid_to_apex(model, &rows, outward0, tolerance, tol)?
    } else {
        let outward0 = outward_at(&rows, &cap_planes, false);
        let outward1 = outward_at(&rows, &cap_planes, true);
        skinned_solid(model, &rows, (outward0, outward1), tolerance, tol)?
    };
    for section in sections {
        built.history.generate(section, built.shape.clone());
    }
    Ok(built)
}

/// Sample a closed wire at `count` matched arc-length fractions.
fn sample_wire(
    model: &Model,
    wire: &Shape,
    count: usize,
    tol: Tolerances,
) -> OgeomResult<Vec<Point>> {
    sample_wire_from(model, wire, count, None, tol)
}

/// As [`sample_wire`], with the arc-length origin rotated to the dense
/// sample nearest `start_hint`: how a caller says which point of each
/// section rows up with which, instead of leaning on traversal starts.
fn sample_wire_from(
    model: &Model,
    wire: &Shape,
    count: usize,
    start_hint: Option<Point>,
    tol: Tolerances,
) -> OgeomResult<Vec<Point>> {
    let mut dense = dense_wire(model, wire, tol)?;
    if let Some(hint) = start_hint {
        let mut best = 0usize;
        let mut held = f64::INFINITY;
        for (i, p) in dense.iter().enumerate() {
            let d = p.distance(hint);
            if d < held {
                held = d;
                best = i;
            }
        }
        dense.rotate_left(best);
    }
    Ok(resampled(&dense, count))
}

/// A closed wire as a dense polyline, by traversal, without the closing
/// point.
fn dense_wire(model: &Model, wire: &Shape, tol: Tolerances) -> OgeomResult<Vec<Point>> {
    let mut dense: Vec<Point> = Vec::new();
    for edge in explore(model, wire, Filter::OfType(ShapeType::Edge))? {
        let Some(data) = model.node(&edge).and_then(|n| n.data().as_edge()) else {
            ogeom_bail!(Construction, "a section edge holds no data");
        };
        let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
            ogeom_bail!(Construction, "a section edge has no curve");
        };
        let Some(geometry) = model.geometry().curve(*curve) else {
            ogeom_bail!(Dangling, "curve is not in this model");
        };
        let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
        // Where the edge stands: a placed wire's edges carry its placement.
        let placement = edge.transform(model.datums())?;
        let at = |f: f64| -> OgeomResult<Point> {
            let t = if reversed {
                range.1 - (range.1 - range.0) * f
            } else {
                range.0 + (range.1 - range.0) * f
            };
            Ok(placement.apply(geometry.point_at(t, tol)?))
        };
        // Samples are read along the chords between these points, so the
        // chords keep within a micron of the curve: the sag measured over a
        // first pass, and the count raised by its square root (the sag
        // falls with the square of the count).
        const FIRST: u32 = 64;
        const MOST: u32 = 8192;
        let mut sag = 0.0_f64;
        for i in 0..FIRST {
            let (f0, f1) = (
                f64::from(i) / f64::from(FIRST),
                f64::from(i + 1) / f64::from(FIRST),
            );
            let (a, b) = (at(f0)?, at(f1)?);
            let mid = at(f64::midpoint(f0, f1))?;
            sag = sag.max(mid.distance(a.midpoint(b)));
        }
        let allowed = tol.confusion() * 10.0;
        #[allow(
            clippy::cast_possible_truncation,
            clippy::cast_sign_loss,
            reason = "a sample count, bounded"
        )]
        let count = if sag <= allowed {
            FIRST
        } else {
            ((f64::from(FIRST) * (sag / allowed).sqrt()).ceil() as u32).clamp(FIRST, MOST)
        };
        for i in 0..count {
            dense.push(at(f64::from(i) / f64::from(count))?);
        }
    }
    if dense.is_empty() {
        ogeom_bail!(Construction, "a section has no edges");
    }
    Ok(dense)
}

/// `count` points at even fractions of a closed polyline's length, from
/// its first point.
fn resampled(dense: &[Point], count: usize) -> Vec<Point> {
    let mut lengths = vec![0.0];
    for w in dense.windows(2) {
        let last = lengths[lengths.len() - 1];
        lengths.push(last + w[0].distance(w[1]));
    }
    let closing = dense[dense.len() - 1].distance(dense[0]);
    let total = lengths[lengths.len() - 1] + closing;
    let mut out = Vec::with_capacity(count);
    let mut cursor = 0usize;
    for s in 0..count {
        #[allow(clippy::cast_precision_loss)]
        let target = total * (s as f64) / (count as f64);
        while cursor + 1 < lengths.len() && lengths[cursor + 1] < target {
            cursor += 1;
        }
        let (a, b) = (dense[cursor], dense[(cursor + 1) % dense.len()]);
        let la = lengths[cursor];
        let lb = if cursor + 1 < lengths.len() {
            lengths[cursor + 1]
        } else {
            total
        };
        let f = if lb > la {
            (target - la) / (lb - la)
        } else {
            0.0
        };
        out.push(a + (b - a) * f.clamp(0.0, 1.0));
    }
    out
}

/// A section's samples in its own frame, started and run where they best
/// match the section before (`previous`, in its frame): a section's start
/// and sense are accidents of how it was drawn, and matched as given they
/// twist the blend. Every start of the dense polyline is tried both ways;
/// the given start and sense win a tie.
fn matched_samples(dense: &[Point], count: usize, previous: &[Point]) -> Vec<Point> {
    let cost = |samples: &[Point]| -> f64 {
        samples
            .iter()
            .zip(previous)
            .map(|(p, q)| (*p - *q).dot(*p - *q))
            .sum()
    };
    let mut best = resampled(dense, count);
    let mut held = cost(&best);
    let scale: f64 = previous
        .iter()
        .chain(&best)
        .map(|p| (*p - Point::ORIGIN).dot(*p - Point::ORIGIN))
        .sum();
    let slack = scale * 1e-12;
    let backward: Vec<Point> = dense.iter().rev().copied().collect();
    for way in [dense, &backward[..]] {
        for start in 0..way.len() {
            let mut turned = way.to_vec();
            turned.rotate_left(start);
            let samples = resampled(&turned, count);
            let c = cost(&samples);
            if c < held - slack {
                held = c;
                best = samples;
            }
        }
    }
    best
}

/// Sweep a circular profile along a free-form spine, skinned.
///
/// Frames along the spine are rotation-minimizing (the double-reflection
/// construction), so the tube neither twists nor kinks where the spine
/// bends; the skin holds the sampled circles to `tolerance`, and the caps
/// sit perpendicular to the spine's ends.
///
/// # Errors
///
/// As [`make_pipe`], plus [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) if
/// the skin cannot reach the
/// tolerance.
pub fn make_pipe_skinned(
    model: &mut Model,
    spine: &Shape,
    radius: f64,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if !radius.is_finite() || radius <= tol.confusion() {
        ogeom_bail!(Construction, "a pipe of radius {radius} holds nothing");
    }
    let (curve, range) = {
        let Some(data) = model.node(spine).and_then(|n| n.data().as_edge()) else {
            ogeom_bail!(Construction, "a pipe runs along an edge");
        };
        let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
            ogeom_bail!(Construction, "the spine has no curve");
        };
        let Some(geometry) = model.geometry().curve(*curve) else {
            ogeom_bail!(Dangling, "curve is not in this model");
        };
        (geometry.clone(), *range)
    };
    const STATIONS: usize = 33;
    const AROUND: usize = 40;
    let mut stations: Vec<SpineStation> = Vec::with_capacity(STATIONS);
    for i in 0..STATIONS {
        #[allow(clippy::cast_precision_loss)]
        let t = range.0 + (range.1 - range.0) * (i as f64) / ((STATIONS - 1) as f64);
        let p = curve.point_at(t, tol)?;
        let d = curve.d1_at(t, tol)?;
        let m = d.magnitude();
        if m <= tol.confusion() {
            ogeom_bail!(Construction, "the spine is degenerate at {t}");
        }
        stations.push(SpineStation {
            at: p,
            tangent: d / m,
            edge: 0,
            t,
        });
    }
    let normals = rmf_normals(&stations);
    let mut rows: Vec<Vec<Point>> = Vec::with_capacity(STATIONS);
    for (i, station) in stations.iter().enumerate() {
        let x = normals[i];
        let y = station.tangent.cross(x);
        let mut row = Vec::with_capacity(AROUND);
        for a in 0..AROUND {
            #[allow(clippy::cast_precision_loss)]
            let ang = core::f64::consts::TAU * (a as f64) / (AROUND as f64);
            row.push(station.at + (x * ang.cos() + y * ang.sin()) * radius);
        }
        rows.push(row);
    }
    let mut built = skinned_solid(
        model,
        &rows,
        (
            EndCap::Plane(-stations[0].tangent),
            EndCap::Plane(stations[STATIONS - 1].tangent),
        ),
        tolerance,
        tol,
    )?;
    built.history.generate(spine, built.shape.clone());
    Ok(built)
}

/// One sampled spine station: where the spine is and which way it runs.
#[derive(Clone, Copy)]
pub(crate) struct SpineStation {
    pub(crate) at: Point,
    /// The unit tangent, in the direction of travel.
    pub(crate) tangent: Vector,
    /// The spine edge this station stands on, by position in the spine.
    pub(crate) edge: usize,
    /// The station's parameter on that edge's curve.
    pub(crate) t: f64,
}

/// One profile wire's closed shell round the spine: smooth wires skin as a
/// single closed face, faceted ones as one ring strip per facet.
#[allow(clippy::too_many_arguments, reason = "one frame, spelled out")]
fn closed_loop_shell(
    model: &mut Model,
    profile_loop: &Shape,
    edges: &[Shape],
    smooth: bool,
    stations: &[SpineStation],
    normals: &[Vector],
    frame0: (Point, Vector),
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    const AROUND: usize = 40;
    let (origin, x0) = frame0;
    let t0 = stations[0].tangent;
    let y0 = t0.cross(x0);
    if !smooth {
        // A faceted profile: one ring strip per profile edge; a fit cannot
        // speak a corner, so each facet gets its own v-closed skin and the
        // strips sew along the corner loops they share within tolerance.
        const ALONG_EDGE: usize = 8;
        // Outward for a ring strip means away from the spine's own line,
        // not from the loop's centroid: a ring's inner side *faces* the
        // centroid. The hint is the station the strip's midpoint rides.
        let mid_station = stations[stations.len() / 2].at;
        let mut faces = Vec::with_capacity(edges.len());
        // Each corner loop is one rail edge shared by the two strips that
        // meet along it, the wrap included.
        let mut rails: Vec<Option<Shape>> = vec![None; edges.len()];
        for (index, edge) in edges.iter().enumerate() {
            let (curve, range) = spine_curve_of(model, edge)?;
            let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
            let mut flat_row: Vec<(f64, f64)> = Vec::with_capacity(ALONG_EDGE + 1);
            for kk in 0..=ALONG_EDGE {
                #[allow(clippy::cast_precision_loss)]
                let f = (kk as f64) / (ALONG_EDGE as f64);
                let t = if reversed {
                    range.1 - (range.1 - range.0) * f
                } else {
                    range.0 + (range.1 - range.0) * f
                };
                let p = curve.point_at(t, tol)?;
                flat_row.push(((p - origin).dot(x0), (p - origin).dot(y0)));
            }
            // rows[j = station][i = across the facet].
            let rows: Vec<Vec<Point>> = stations
                .iter()
                .enumerate()
                .map(|(i, station)| {
                    let x = normals[i];
                    let y = station.tangent.cross(x);
                    flat_row
                        .iter()
                        .map(|(a, b)| station.at + x * *a + y * *b)
                        .collect()
                })
                .collect();
            let next = (index + 1) % edges.len();
            let shared = [rails[index].clone(), rails[next].clone()];
            let (face, rail0, rail1) = skinned_ring_strip(
                model,
                &rows,
                mid_station,
                [shared[0].as_ref(), shared[1].as_ref()],
                tolerance,
                tol,
            )?;
            rails[index] = Some(rail0);
            rails[next] = Some(rail1);
            faces.push(face);
        }
        let sewn = sew(model, &faces, tol)?;
        if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
            if std::env::var_os("OGEOM_DEBUG_RING").is_some() {
                use ogeom_geom::Curve3d as _;
                eprintln!(
                    "RING: {} shells from {} strips",
                    sewn.shells.len(),
                    faces.len()
                );
                for shell in &sewn.shells {
                    for edge in ogeom_topo::explore_unique(model, shell, ShapeType::Edge)? {
                        let mut uses = 0;
                        for f in explore(model, shell, Filter::OfType(ShapeType::Face))? {
                            for w in model.children_of(&f)? {
                                for e in model.children_of(&w)? {
                                    if e.node() == edge.node() {
                                        uses += 1;
                                    }
                                }
                            }
                        }
                        if uses == 1
                            && let Some(d) = model.node(&edge).and_then(|n| n.data().as_edge())
                            && let Some(ogeom_topo::EdgeRepr::Curve3d { curve, range, .. }) =
                                d.curve3d()
                            && let Some(g) = model.geometry().curve(*curve)
                        {
                            eprintln!(
                                "RING open edge tol {:.2e}: {:?} -> {:?}",
                                d.tolerance.get(),
                                g.point_at(range.0, tol)?,
                                g.point_at(range.1, tol)?
                            );
                        }
                    }
                }
            }
            ogeom_bail!(Construction, "the faceted ring did not close");
        }
        return Ok(sewn.shells[0].clone());
    }
    let samples = sample_wire(model, profile_loop, AROUND, tol)?;
    let flat: Vec<(f64, f64)> = samples
        .iter()
        .map(|p| ((*p - origin).dot(x0), (*p - origin).dot(y0)))
        .collect();
    let rows: Vec<Vec<Point>> = stations
        .iter()
        .enumerate()
        .map(|(i, station)| {
            let x = normals[i];
            let y = station.tangent.cross(x);
            flat.iter()
                .map(|(a, b)| station.at + x * *a + y * *b)
                .collect()
        })
        .collect();
    closed_skinned_shell(model, &rows, tolerance, tol)
}

/// Rotation-minimizing normals along the stations, by double reflection:
/// reflect in each chord's plane, then in the plane bisecting the tangents.
/// Self-contained (it needs only the station list) and shared by every
/// sweep that must not twist where its spine bends.
fn rmf_normals(stations: &[SpineStation]) -> Vec<Vector> {
    let mut normals: Vec<Vector> = Vec::with_capacity(stations.len());
    let t0 = stations[0].tangent;
    let seed = if t0.cross(ogeom_math::Vector::Z).magnitude() > 0.5 {
        ogeom_math::Vector::Z
    } else {
        ogeom_math::Vector::X
    };
    let n0 = {
        let v = seed - t0 * seed.dot(t0);
        v / v.magnitude()
    };
    normals.push(n0);
    for i in 1..stations.len() {
        let (p0, t0) = (stations[i - 1].at, stations[i - 1].tangent);
        let (p1, t1) = (stations[i].at, stations[i].tangent);
        let n = normals[i - 1];
        let v1 = p1 - p0;
        let c1 = v1.dot(v1);
        if c1 <= 1e-20 {
            // A corner's twin station: no travel to reflect through. The
            // normal is reflected across the corner's mitre plane instead;
            // for a vector square to the incoming tangent that is exactly
            // the parallel transport about the corner's own axis, and the
            // mirror symmetry is what lands both legs' sheared sections on
            // one ring. A planar corner's normal lies in the mitre plane
            // already and carries straight across; a skew corner's does
            // not, and carried unchanged it leaves the far leg's section
            // off the mitre.
            let bisector = t0 + t1;
            let m = bisector.magnitude();
            if m <= 1e-12 {
                normals.push(n);
                continue;
            }
            let b = bisector / m;
            normals.push(n - b * (2.0 * n.dot(b)));
            continue;
        }
        normals.push(rmf_step(p0, t0, n, p1, t1));
    }
    normals
}

/// One rotation-minimizing step: the normal `n0` at `(p0, t0)` carried to
/// `(p1, t1)` by double reflection. No travel means no change.
fn rmf_step(p0: Point, t0: Vector, n0: Vector, p1: Point, t1: Vector) -> Vector {
    let v1 = p1 - p0;
    let c1 = v1.dot(v1);
    if c1 <= 1e-20 {
        return n0;
    }
    let nl = n0 - v1 * (2.0 / c1 * v1.dot(n0));
    let tl = t0 - v1 * (2.0 / c1 * v1.dot(t0));
    let v2 = t1 - tl;
    let c2 = v2.dot(v2);
    let next = if c2 > 1e-20 {
        nl - v2 * (2.0 / c2 * v2.dot(nl))
    } else {
        nl
    };
    next / next.magnitude()
}

/// A leg's generators, evaluated anywhere: the spine's own curve between
/// stations with the rotation-minimizing normal carried one step from the
/// station behind, and a straight extension past either end in the end
/// frame: the surface a mitre trims against. Parameters are station
/// indices; a unit beyond an end is one station spacing.
struct SpineWalk<'a> {
    /// Each spine edge's curve, range and whether it is travelled reversed.
    curves: Vec<(ogeom_geom::Curve, (f64, f64), bool)>,
    stations: &'a [SpineStation],
    normals: &'a [Vector],
}

/// Where two legs' generators for one profile point meet at a corner: the
/// point, each leg's parameter, and how far the two generators actually
/// miss each other (zero when the corner turns in the plane).
struct CornerJoin {
    at: Point,
    s1: f64,
    s2: f64,
    gap: f64,
}

impl SpineWalk<'_> {
    /// The spine point, unit tangent and frame normal at `s` within the run
    /// `(rs, re)`.
    fn frame_at(
        &self,
        s: f64,
        (rs, re): (usize, usize),
        tol: Tolerances,
    ) -> OgeomResult<(Point, Vector, Vector)> {
        let st = self.stations;
        let at = |i: usize| (st[i].at, st[i].tangent, self.normals[i]);
        #[allow(clippy::cast_precision_loss)]
        let (rsf, ref_) = (rs as f64, re as f64);
        if s <= rsf {
            let (p, t, n) = at(rs);
            let h = st[rs].at.distance(st[(rs + 1).min(re)].at);
            return Ok((p + t * ((s - rsf) * h), t, n));
        }
        if s >= ref_ {
            let (p, t, n) = at(re);
            let h = st[re].at.distance(st[re.saturating_sub(1).max(rs)].at);
            return Ok((p + t * ((s - ref_) * h), t, n));
        }
        #[allow(clippy::cast_possible_truncation, clippy::cast_sign_loss)]
        let j = s.floor() as usize;
        #[allow(clippy::cast_precision_loss)]
        let f = s - j as f64;
        if f <= 0.0 {
            return Ok(at(j));
        }
        let (curve, range, reversed) = &self.curves[st[j + 1].edge];
        let t_from = if st[j].edge == st[j + 1].edge {
            st[j].t
        } else if *reversed {
            range.1
        } else {
            range.0
        };
        let t = t_from + (st[j + 1].t - t_from) * f;
        let p = curve.point_at(t, tol)?;
        let d = curve.d1_at(t, tol)?;
        let m = d.magnitude();
        if m <= tol.confusion() {
            ogeom_bail!(Construction, "the spine is degenerate at {t}");
        }
        let tangent = if *reversed { -(d / m) } else { d / m };
        let n = rmf_step(st[j].at, st[j].tangent, self.normals[j], p, tangent);
        Ok((p, tangent, n))
    }

    /// The generator of profile point `(a, b)` at `s` within the run.
    fn generator(
        &self,
        s: f64,
        run: (usize, usize),
        (a, b): (f64, f64),
        tol: Tolerances,
    ) -> OgeomResult<Point> {
        let (p, t, x) = self.frame_at(s, run, tol)?;
        let y = t.cross(x);
        Ok(p + x * a + y * b)
    }

    /// Where the generators of one profile point on the leg `before` and
    /// the leg `after` a corner meet: Gauss-Newton on both parameters from
    /// the corner itself, minimising the distance between the two.
    fn join(
        &self,
        before: (usize, usize),
        after: (usize, usize),
        ab: (f64, f64),
        tol: Tolerances,
    ) -> OgeomResult<CornerJoin> {
        const STEP: f64 = 1e-4;
        #[allow(clippy::cast_precision_loss)]
        let (mut s1, mut s2) = (before.1 as f64, after.0 as f64);
        for _ in 0..60 {
            let g1 = self.generator(s1, before, ab, tol)?;
            let g2 = self.generator(s2, after, ab, tol)?;
            let f = g1 - g2;
            let d1 = (self.generator(s1 + STEP, before, ab, tol)?
                - self.generator(s1 - STEP, before, ab, tol)?)
                / (2.0 * STEP);
            let d2 = (self.generator(s2 + STEP, after, ab, tol)?
                - self.generator(s2 - STEP, after, ab, tol)?)
                / (2.0 * STEP);
            let (a11, a12, a22) = (d1.dot(d1), -d1.dot(d2), d2.dot(d2));
            let (b1, b2) = (-f.dot(d1), f.dot(d2));
            let det = a11 * a22 - a12 * a12;
            if det.abs() <= 1e-30 {
                break;
            }
            let e1 = (b1 * a22 - a12 * b2) / det;
            let e2 = (a11 * b2 - a12 * b1) / det;
            s1 += e1;
            s2 += e2;
            if e1.abs().max(e2.abs()) <= 1e-12 {
                break;
            }
        }
        let g1 = self.generator(s1, before, ab, tol)?;
        let g2 = self.generator(s2, after, ab, tol)?;
        Ok(CornerJoin {
            at: g1.midpoint(g2),
            s1,
            s2,
            gap: g1.distance(g2),
        })
    }
}

/// Sweep a planar profile along a helix about `axis`: a screw motion,
/// every point of the profile running its own helix. The thread and spring
/// operation, with the profile in a plane through the axis; a profile in
/// any other plane clear of the axis (square to it, a ramp or a stair's
/// tread) climbs the same way.
///
/// `pitch` is the advance per turn along `axis`, `turns` how far the
/// profile turns, `left_handed` turns it the other way about the axis for
/// the same advance, and `taper_per_turn` moves every point away from the
/// axis by that much per turn (a conical helix; zero for a cylindrical
/// one). A pitch of zero with a taper is a flat spiral, every point
/// turning in its plane square to the axis while moving out. The walls are
/// fitted through each profile edge's exact screw images; the caps are the
/// profile where it starts and where it ends.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if
/// the profile is not a planar face, meets the axis, would meet itself one
/// turn on (its extent along the axis is not less than the pitch, or for a
/// flat spiral its extent away from the axis not less than the taper), is
/// carried partly forward through its own plane and partly back, or tapers
/// onto the axis; if a taper is asked of a profile whose plane neither
/// holds the axis nor stands square to it; if `turns` is not positive,
/// `pitch` is negative, or both `pitch` and the taper are zero.
/// [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) if a wall
/// cannot be fitted.
#[allow(clippy::too_many_arguments, clippy::too_many_lines)]
pub fn make_helical_sweep(
    model: &mut Model,
    profile: &Shape,
    axis: ogeom_math::Axis,
    pitch: f64,
    turns: f64,
    left_handed: bool,
    taper_per_turn: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if !(pitch.is_finite() && pitch >= 0.0) || !(turns.is_finite() && turns > 0.0) {
        ogeom_bail!(
            Construction,
            "a helical sweep needs a pitch of zero or more and a positive turn count; \
             got {pitch} and {turns}"
        );
    }
    if !taper_per_turn.is_finite() {
        ogeom_bail!(Construction, "a taper of {taper_per_turn} is not a length");
    }
    let flat = pitch <= tol.confusion();
    if flat && taper_per_turn.abs() <= tol.confusion() {
        ogeom_bail!(
            Construction,
            "a helical sweep with no pitch and no taper turns the profile onto itself"
        );
    }
    if model.kind_of(profile)? != ShapeType::Face {
        ogeom_bail!(Construction, "a helical sweep sweeps a planar face");
    }
    let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
        ogeom_bail!(Construction, "a helical sweep sweeps a planar face");
    };
    let z = axis.direction.vector();
    let n = plane.normal().vector();
    let holds_axis = n.dot(z).abs() <= tol.angular()
        && plane.distance_to(axis.location) <= tol.confusion() * 100.0;
    let level = (n.dot(z).abs() - 1.0).abs() <= tol.angular();
    // A taper moves every point away from the axis by the same amount, which
    // keeps a profile in a plane only where that plane holds the axis or
    // stands square to it.
    if taper_per_turn.abs() > tol.confusion() && !holds_axis && !level {
        ogeom_bail!(
            Construction,
            "a tapered helical sweep carries a profile in a plane through the \
             axis or square to it; an oblique profile would leave its plane"
        );
    }
    // A profile the axis runs through has points with no helix to follow:
    // where the axis pierces the plane is asked of the profile's rings,
    // read in the plane.
    if !holds_axis && n.dot(z).abs() > tol.angular() {
        let reach = (plane.origin() - axis.location).dot(n) / z.dot(n);
        let pierce = axis.location + z * reach;
        let frame = plane.frame();
        let flat = |p: Point| {
            let l = frame.to_local(p);
            Point2::new(l.x, l.y)
        };
        let at = flat(pierce);
        let mut inside = false;
        let mut touches = false;
        for wire in explore(model, profile, Filter::OfType(ShapeType::Wire))? {
            let ring: Vec<Point2> = sample_wire(model, &wire, 256, tol)?
                .into_iter()
                .map(flat)
                .collect();
            for (k, a) in ring.iter().enumerate() {
                let b = ring[(k + 1) % ring.len()];
                let (ab, ap) = (b - *a, at - *a);
                let t = (ap.dot(ab) / ab.dot(ab).max(f64::MIN_POSITIVE)).clamp(0.0, 1.0);
                touches |= (*a + ab * t).distance(at) <= tol.confusion() * 100.0;
                if (a.y > at.y) != (b.y > at.y) && at.x < a.x + (at.y - a.y) / (b.y - a.y) * ab.x {
                    inside = !inside;
                }
            }
        }
        if inside || touches {
            ogeom_bail!(
                Construction,
                "the profile meets the axis, where a helical sweep has no helix \
                 to follow"
            );
        }
    }
    let total = core::f64::consts::TAU * turns;
    let sense = if left_handed { -1.0 } else { 1.0 };
    // How a point sets off at the start of the turn: up the axis by the
    // pitch, round it, and out by the taper, each per radian.
    let travel = |p: Point| -> Vector {
        let foot = axis.project(p);
        let out = p - foot;
        let rho = out.magnitude().max(f64::MIN_POSITIVE);
        let radial = out / rho;
        (z * pitch + radial * taper_per_turn) / core::f64::consts::TAU
            + z.cross(radial) * (rho * sense)
    };
    // The screw image of a point after turning through `theta`.
    let screw = |p: Point, theta: f64| -> OgeomResult<Point> {
        let foot = axis.project(p);
        let out = p - foot;
        let rho = out.magnitude();
        let grown = rho + taper_per_turn * theta / core::f64::consts::TAU;
        if rho <= tol.confusion() || grown <= tol.confusion() {
            ogeom_bail!(
                Construction,
                "the profile reaches the axis, where a helical sweep has no \
                 helix to follow"
            );
        }
        let radial = out / rho;
        let across = z.cross(radial);
        let (sin, cos) = (sense * theta).sin_cos();
        let turned = radial * cos + across * sin;
        Ok(foot + z * (pitch * theta / core::f64::consts::TAU) + turned * grown)
    };

    // The profile's loops, their edges in ring order, each sampled.
    let loops = explore(model, profile, Filter::OfType(ShapeType::Wire))?;
    if loops.is_empty() {
        ogeom_bail!(Construction, "the profile has no loop to sweep");
    }
    // Every point leaves the profile's plane the same way, or the profile
    // sweeps back through where it has been; a plane the motion runs
    // along everywhere (a level profile with no pitch) sweeps nothing.
    {
        let (mut ahead, mut behind) = (0.0_f64, 0.0_f64);
        for wire in &loops {
            for p in sample_wire(model, wire, 64, tol)? {
                let t = travel(p);
                let across = t.dot(n) / t.magnitude().max(f64::MIN_POSITIVE);
                ahead = ahead.max(across);
                behind = behind.max(-across);
            }
        }
        if ahead > tol.angular() && behind > tol.angular() {
            ogeom_bail!(
                Construction,
                "the screw carries part of the profile forward through its plane \
                 and part back; the sweep runs into itself"
            );
        }
        if ahead <= tol.angular() && behind <= tol.angular() {
            ogeom_bail!(
                Construction,
                "the screw carries the profile along its own plane; it sweeps \
                 no volume"
            );
        }
    }
    // One turn on, the profile must clear itself: along the axis by the
    // pitch, or for a flat spiral away from it by the taper.
    if turns > 1.0 {
        let mut low = f64::INFINITY;
        let mut high = f64::NEG_INFINITY;
        for wire in &loops {
            for p in sample_wire(model, wire, 64, tol)? {
                let h = if flat {
                    p.distance(axis.project(p))
                } else {
                    (p - axis.location).dot(z)
                };
                low = low.min(h);
                high = high.max(h);
            }
        }
        let (clearance, across) = if flat {
            (taper_per_turn.abs(), "away from the axis")
        } else {
            (pitch, "along the axis")
        };
        if high - low >= clearance - tol.confusion() {
            ogeom_bail!(
                Construction,
                "the profile spans {} {across}, not less than the turn's advance \
                 {clearance}; a turn on it meets itself",
                high - low
            );
        }
    }
    let tolerance = tol.confusion() * 100.0;
    // The walls in quarter turns, each a strip of its own sharing its
    // borders with the next: one fit down many turns of a helix cannot
    // reach the tolerance, a quarter turn's can. The borders stand an
    // eighth of a turn off the profile's own plane and the planes square to
    // it: a strip's border is a copy of the profile, and a plane through the
    // axis at those angles (a block's face set on the axis) would otherwise
    // meet a strip only along its border, where no section starts.
    let quarter = core::f64::consts::FRAC_PI_2;
    let mut borders = vec![0.0];
    let mut next_border = quarter / 2.0;
    while next_border < total - quarter * 1e-3 {
        borders.push(next_border);
        next_border += quarter;
    }
    borders.push(total);
    let segments = borders.len() - 1;
    // The fit keeps fewer controls than samples, so its reach at the
    // samples is set by how many there are.
    const PER_SEGMENT: usize = 48;
    #[allow(clippy::cast_precision_loss)]
    let theta_at = |seg: usize, i: usize| {
        borders[seg] + (borders[seg + 1] - borders[seg]) * (i as f64) / (PER_SEGMENT as f64)
    };

    let mut faces: Vec<Shape> = Vec::new();
    let mut cap_loops: [Vec<Vec<Shape>>; 2] = [Vec::new(), Vec::new()];
    for (li, wire) in loops.iter().enumerate() {
        let hole = li != 0;
        let edges = model.ordered_children_of(wire)?;
        let centre = {
            let samples = sample_wire(model, wire, 32, tol)?;
            #[allow(clippy::cast_precision_loss)]
            let n = samples.len() as f64;
            let sum = samples
                .iter()
                .fold(Vector::new(0.0, 0.0, 0.0), |acc, p| acc + p.to_vector());
            Point::from_vector(sum / n)
        };
        // The ring's pieces: every edge in the ring's sense, a closed one
        // (a circle, the whole ring) cut in quarters so each strip's fit
        // spans a quarter turn round it at most.
        let mut pieces: Vec<(ogeom_geom::Curve, f64, f64)> = Vec::new();
        for edge in &edges {
            let (curve, range) = spine_curve_of(model, edge)?;
            let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
            let (t0, t1) = if reversed { (range.1, range.0) } else { range };
            let closed =
                ogeom_algo::edge_vertices(model, edge)?.is_some_and(|(a, b)| a.is_same(&b));
            let parts = if closed { 4 } else { 1 };
            for k in 0..parts {
                #[allow(clippy::cast_precision_loss)]
                let (f0, f1) = (k as f64 / parts as f64, (k + 1) as f64 / parts as f64);
                pieces.push((curve.clone(), t0 + (t1 - t0) * f0, t0 + (t1 - t0) * f1));
            }
        }
        // Each piece's samples, in the ring's sense, ends on its corners
        // exactly.
        let starts: Vec<Point> = pieces
            .iter()
            .map(|(curve, a, _)| curve.point_at(*a, tol))
            .collect::<OgeomResult<_>>()?;
        let count = pieces.len();
        let mut rows0: Vec<Vec<Point>> = Vec::with_capacity(count);
        for (pi, (curve, a, b)) in pieces.iter().enumerate() {
            let along = if matches!(curve, ogeom_geom::Curve::Line(_)) {
                8
            } else {
                24
            };
            let mut row = Vec::with_capacity(along + 1);
            for k in 0..=along {
                #[allow(clippy::cast_precision_loss)]
                let f = (k as f64) / (along as f64);
                row.push(curve.point_at(a + (b - a) * f, tol)?);
            }
            row[0] = starts[pi];
            row[along] = starts[(pi + 1) % count];
            rows0.push(row);
        }
        let rows_of = |row0: &[Point], seg: usize| -> OgeomResult<Vec<Vec<Point>>> {
            (0..=PER_SEGMENT)
                .map(|i| {
                    row0.iter()
                        .map(|p| screw(*p, theta_at(seg, i)))
                        .collect::<OgeomResult<Vec<Point>>>()
                })
                .collect()
        };

        // A vertex set at every segment boundary.
        let mut corners: Vec<Vec<Shape>> = Vec::with_capacity(segments + 1);
        for b in 0..=segments {
            let theta = if b == segments { total } else { theta_at(b, 0) };
            let mut set = Vec::with_capacity(count);
            for p in &starts {
                set.push(ogeom_algo::make_vertex(model, screw(*p, theta)?).shape);
            }
            corners.push(set);
        }
        let mut bottoms = Vec::with_capacity(count);
        let mut tops = Vec::with_capacity(count);
        // Per edge, the previous segment's top border, for the next to
        // start on.
        let mut held_tops: Vec<Option<Shape>> = vec![None; count];
        for seg in 0..segments {
            let hint = screw(centre, theta_at(seg, PER_SEGMENT / 2))?;
            let mut first_rail: Option<Shape> = None;
            let mut prev_rail: Option<Shape> = None;
            for ei in 0..count {
                let rows = rows_of(&rows0[ei], seg)?;
                let next = (ei + 1) % count;
                let last_rail = if ei + 1 == count {
                    first_rail.clone()
                } else {
                    None
                };
                let (from, to) = (&corners[seg], &corners[seg + 1]);
                let strip = skinned_strip(
                    model,
                    &rows,
                    (&from[ei], &from[next], &to[ei], &to[next]),
                    [
                        held_tops[ei].as_ref(),
                        None,
                        prev_rail.as_ref(),
                        last_rail.as_ref(),
                    ],
                    hint,
                    hole,
                    false,
                    tolerance,
                    tol,
                )?;
                if ei == 0 {
                    first_rail = Some(strip.rail0.clone());
                }
                prev_rail = Some(strip.rail1.clone());
                faces.push(strip.face.clone());
                if seg == 0 {
                    bottoms.push(strip.bottom.clone());
                }
                if seg + 1 == segments {
                    tops.push(strip.top.clone());
                }
                held_tops[ei] = Some(strip.top);
            }
        }
        cap_loops[0].push(bottoms);
        cap_loops[1].push(tops);
    }

    // The caps: the profile's plane where it starts, and that plane
    // screwed on to where it ends.
    for (end, loops) in cap_loops.iter().enumerate() {
        let theta = if end == 0 { 0.0 } else { total };
        let centre = centre_of(model, profile, tol)?;
        let at = screw(centre, theta)?;
        // The profile's plane, turned with it: out of the solid, back
        // against the motion at the start and on with it at the end.
        let normal = {
            let (sin, cos) = (sense * theta).sin_cos();
            let along = z * n.dot(z);
            let square = n - along;
            let turned = along + square * cos + z.cross(square) * sin;
            let forward = if turned.dot(travel(centre)) >= 0.0 {
                turned
            } else {
                -turned
            };
            if end == 0 { -forward } else { forward }
        };
        let cap_plane = Plane::through(at, Direction::new(normal, tol)?);
        let mut reach = 1.0_f64;
        for edges in loops {
            for edge in edges {
                let (curve, range) = spine_curve_of(model, edge)?;
                for k in 0..8 {
                    let p =
                        curve.point_at(range.0 + (range.1 - range.0) * f64::from(k) / 8.0, tol)?;
                    reach = reach.max(p.distance(at) * 2.0);
                }
            }
        }
        let surface: SurfaceGeometry =
            PlaneSurface::over(cap_plane, (-reach, reach), (-reach, reach))?.into();
        let mut wires = Vec::with_capacity(loops.len());
        for edges in loops {
            wires.push(ogeom_algo::make_wire(model, edges, tol)?.shape);
        }
        let face = ogeom_algo::make_face(model, surface, &wires, tol)?.shape;
        let cap_id = {
            let Some(ogeom_topo::NodeData::Face(data)) = model.node(&face).map(|n| n.data()) else {
                ogeom_bail!(Construction, "the cap holds no face data");
            };
            data.surface
        };
        let frame = cap_plane.frame();
        for edges in loops {
            for edge in edges {
                let (curve, range) = spine_curve_of(model, edge)?;
                let ogeom_geom::Curve::BSpline(bs) = &curve else {
                    ogeom_bail!(Construction, "a swept ring is not a spline");
                };
                let control2: Vec<Point2> = bs
                    .control_points()
                    .iter()
                    .map(|w| {
                        let local = frame.to_local(w.point());
                        Point2::new(local.x, local.y)
                    })
                    .collect();
                let pcurve: ogeom_geom::PlanarCurve =
                    ogeom_geom::BSpline2d::new(bs.knots().clone(), control2, tol)?.into();
                ogeom_algo::attach_pcurve(
                    model,
                    edge,
                    pcurve,
                    cap_id,
                    ogeom_topo::Location::identity(),
                    range,
                )?;
            }
        }
        faces.push(face);
    }

    let sewn = sew(model, &faces, tol)?;
    if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
        ogeom_bail!(Construction, "the helical sweep did not close");
    }
    let solid = make_solid(model, &sewn.shells)?.shape;
    let mut history = History::new();
    history.generate(profile, solid.clone());
    Ok(Built::new(solid, history))
}

/// Revolve `profile` (a planar face) about `axis`, each point turning the
/// right-handed way about it until its circle first meets the surface of
/// `limit`, a face taken as its whole surface the way
/// [`make_half_space`](ogeom_algo::make_half_space) takes it: the revolution
/// "up to face", stopping on the target at a different angle for each
/// point.
///
/// Built as the half turn of the profile kept on the profile's side of the
/// limit: exact wherever every circle meets the surface within that half
/// turn and does not come back to the profile's side before it ends.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the
/// profile is not a planar face whose plane holds the axis, a point's
/// circle never meets the surface, or one meets it only past half a turn or
/// crosses it back within the half turn; as the revolution, the half space
/// and the boolean otherwise.
pub fn make_revolution_until(
    model: &mut Model,
    profile: &Shape,
    axis: ogeom_math::Axis,
    limit: &Shape,
    tol: Tolerances,
) -> OgeomResult<Built> {
    use ogeom_geom::Surface as _;
    if model.kind_of(profile)? != ShapeType::Face {
        ogeom_bail!(
            Construction,
            "a revolution up to a face revolves a planar face"
        );
    }
    if model.kind_of(limit)? != ShapeType::Face {
        ogeom_bail!(Construction, "a revolution stops on a face");
    }
    let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
        ogeom_bail!(
            Construction,
            "a revolution up to a face revolves a planar face"
        );
    };
    let z = axis.direction.vector();
    if plane.normal().vector().dot(z).abs() > tol.angular()
        || plane.distance_to(axis.location) > tol.confusion() * 100.0
    {
        ogeom_bail!(
            Construction,
            "the profile's plane does not hold the axis; a revolution up to a \
             face turns a profile about an axis in its own plane"
        );
    }
    // The limit's surface, placed, and which side of it a point stands.
    let (surface, flip) = {
        let Some(data) = model.node(limit).and_then(|n| n.data().as_face()) else {
            ogeom_bail!(Construction, "the limit holds no face data");
        };
        let Some(surface) = model.geometry().surface(data.surface) else {
            ogeom_bail!(Dangling, "the limit's surface is not in this model");
        };
        (
            surface
                .clone()
                .transformed(&limit.transform(model.datums())?, tol)?,
            limit.orientation() == ogeom_topo::Orientation::Reversed,
        )
    };
    let side = |p: Point| -> OgeomResult<f64> {
        let foot = ogeom_algo::project_on_surface(&surface, p, 16, tol)?;
        let (u, v) = foot.parameters;
        let n = surface.normal_at(u, v, tol)?.vector();
        let s = (p - foot.point).dot(n);
        Ok(if flip { -s } else { s })
    };
    // Every boundary point's circle: its first meeting within half a turn,
    // and no return to the profile's side before the half turn ends.
    let turn_of = |p: Point, theta: f64| Transform::rotation(axis, theta).apply(p);
    let mut rings = Vec::new();
    for wire in explore(model, profile, Filter::OfType(ShapeType::Wire))? {
        rings.extend(sample_wire(model, &wire, 64, tol)?);
    }
    let start_side = side(centre_of(model, profile, tol)?)?.signum();
    const STEPS: u32 = 360;
    for p in &rings {
        let mut met: Option<f64> = None;
        for k in 1..=STEPS {
            let theta = core::f64::consts::TAU * f64::from(k) / f64::from(STEPS);
            let here = side(turn_of(*p, theta))?;
            match met {
                None if here.signum() != start_side && here.abs() > tol.confusion() => {
                    met = Some(theta);
                }
                Some(first)
                    if theta <= core::f64::consts::PI
                        && here.signum() == start_side
                        && here.abs() > tol.confusion() =>
                {
                    ogeom_bail!(
                        Construction,
                        "the circle through {p:?} meets the limit at {first} and \
                         crosses back within half a turn; a revolution up to it \
                         is not built there"
                    );
                }
                _ => {}
            }
        }
        match met {
            None => ogeom_bail!(
                Construction,
                "the circle through {p:?} never meets the limit's surface"
            ),
            Some(first) if first > core::f64::consts::PI => ogeom_bail!(
                Construction,
                "the circle through {p:?} meets the limit only past half a turn"
            ),
            _ => {}
        }
    }
    let half = ogeom_algo::make_revolution(model, profile, axis, core::f64::consts::PI, tol)?;
    let inside = centre_of(model, profile, tol)?;
    let bound = ogeom_algo::make_half_space(model, limit, inside, tol)?.shape;
    let mut built = ogeom_bool::common(model, &half.shape, &bound, tol)?;
    built.history.generate(profile, built.shape.clone());
    Ok(built)
}

/// The centroid of a face's outer ring's samples.
fn centre_of(model: &Model, profile: &Shape, tol: Tolerances) -> OgeomResult<Point> {
    let Some(wire) = explore(model, profile, Filter::OfType(ShapeType::Wire))?
        .into_iter()
        .next()
    else {
        ogeom_bail!(Construction, "the profile has no loop");
    };
    let samples = sample_wire(model, &wire, 32, tol)?;
    #[allow(clippy::cast_precision_loss)]
    let n = samples.len() as f64;
    let sum = samples
        .iter()
        .fold(Vector::new(0.0, 0.0, 0.0), |acc, p| acc + p.to_vector());
    Ok(Point::from_vector(sum / n))
}

/// The pipe along a spine of lines and circular arcs meeting tangent to
/// one another, built exactly: down a line the section is extruded, round
/// an arc it is revolved about the arc's axis (the rotation-minimizing
/// frame of a circle is its own rotation), and the legs are fused on the
/// sections they share. Every wall is then the closed form its profile
/// edge sweeps: a plane, drum, cone, ball or torus where the edge is a
/// line or circle. `None` where the spine or profile is not of that kind,
/// for the general construction to take.
fn exact_legs(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    frenet: bool,
    tol: Tolerances,
) -> OgeomResult<Option<Built>> {
    use ogeom_geom::{Curve, Curve3d as _};
    let edges: Vec<Shape> = match model.kind_of(spine)? {
        ShapeType::Edge => vec![spine.clone()],
        ShapeType::Wire => model.ordered_children_of(spine)?,
        _ => return Ok(None),
    };
    let solid = match model.kind_of(profile)? {
        ShapeType::Face => true,
        ShapeType::Wire => false,
        _ => return Ok(None),
    };
    if edges.is_empty() || (!solid && edges.len() > 1) {
        return Ok(None);
    }
    // Each leg: where it starts and ends, its heading at either end, and
    // the motion that carries the section down it.
    struct Leg {
        start: Point,
        heading: (Vector, Vector),
        motion: Transform,
        along: LegKind,
    }
    enum LegKind {
        Line(Vector),
        Arc(ogeom_math::Axis, f64),
    }
    let mut legs = Vec::with_capacity(edges.len());
    for edge in &edges {
        let (curve, range) = spine_curve_of(model, edge)?;
        let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
        let (t0, t1) = if reversed { (range.1, range.0) } else { range };
        let (a, b) = (curve.point_at(t0, tol)?, curve.point_at(t1, tol)?);
        let sense = if reversed { -1.0 } else { 1.0 };
        let heading = (curve.d1_at(t0, tol)? * sense, curve.d1_at(t1, tol)? * sense);
        let basis = match &curve {
            Curve::Trimmed(t) => t.basis().clone(),
            other => other.clone(),
        };
        let (motion, along) = match basis {
            Curve::Line(_) => {
                if frenet {
                    return Ok(None);
                }
                (Transform::translation(b - a), LegKind::Line(b - a))
            }
            Curve::Circle(c) => {
                let frame = c.circle().frame();
                // Turning the way the walk heads at its start.
                let turn = (a - frame.origin()).cross(heading.0);
                let direction = if turn.dot(frame.z().vector()) > 0.0 {
                    frame.z()
                } else {
                    -frame.z()
                };
                let axis = ogeom_math::Axis {
                    location: frame.origin(),
                    direction,
                };
                let angle = (t1 - t0).abs();
                (Transform::rotation(axis, angle), LegKind::Arc(axis, angle))
            }
            _ => return Ok(None),
        };
        legs.push(Leg {
            start: a,
            heading,
            motion,
            along,
        });
    }
    // Legs meet tangent, or all are straight and their corners are mitred
    // exactly below; any other corner is mitred by the general path.
    let cornered = legs.windows(2).any(|pair| {
        let (x, y) = (pair[0].heading.1, pair[1].heading.0);
        x.cross(y).magnitude() > tol.angular() * x.magnitude() * y.magnitude() || x.dot(y) <= 0.0
    });
    if cornered {
        let lines: Option<Vec<(Point, Vector)>> = legs
            .iter()
            .map(|leg| match leg.along {
                LegKind::Line(v) => Some((leg.start, v)),
                LegKind::Arc(..) => None,
            })
            .collect();
        return match lines {
            Some(lines) if solid && !frenet => mitred_lines(model, profile, spine, &lines, tol),
            _ => Ok(None),
        };
    }
    // The profile square to the spine's exact start tangent, and on it.
    let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
        return Ok(None);
    };
    let t0 = legs[0].heading.0;
    if plane.normal().vector().cross(t0).magnitude() > tol.angular() * t0.magnitude()
        || plane.distance_to(legs[0].start) > tol.confusion() * 100.0
    {
        return Ok(None);
    }

    let mut joined: Vec<Shape> = Vec::new();
    let mut carried = Transform::IDENTITY;
    for leg in &legs {
        // Rebuilt even where it stands: the caps are this face, and it
        // carries its trims on its own plane.
        let heading = carried.apply_vector(legs[0].heading.0);
        let section = realized_profile_wound(model, profile, &carried, Some(heading), tol)?;
        let piece = match leg.along {
            LegKind::Line(v) => ogeom_algo::make_prism(model, &section, v, tol)?.shape,
            LegKind::Arc(axis, angle) => {
                ogeom_algo::make_revolution(model, &section, axis, angle, tol)?.shape
            }
        };
        carried = leg.motion * carried;
        joined.push(piece);
    }
    let result = fuse_in_order(model, joined, tol)?;
    Ok(result.map(|shape| {
        let mut history = History::new();
        history.generate(spine, shape.clone());
        history.generate(profile, shape.clone());
        Built::new(shape, history)
    }))
}

/// A face swept down a spine of straight legs with corners, exactly: each
/// leg is the profile's prism, run on past its corners and trimmed by the
/// mitre plane that halves each corner, and the pieces are fused. The
/// section turns at each corner by the least rotation taking one leg's
/// heading to the next's, as the rotation-minimizing frame does, and every
/// wall is the closed form its profile edge sweeps along a line.
///
/// `None` where the spine closes on itself (a ring's section may come back
/// turned, which the general construction settles), a corner all but
/// doubles back, or a leg is too short for the corners at its ends to trim
/// it apart.
fn mitred_lines(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    legs: &[(Point, Vector)],
    tol: Tolerances,
) -> OgeomResult<Option<Built>> {
    if let (Some(first), Some(last)) = (legs.first(), legs.last())
        && first.0.distance(last.0 + last.1) <= tol.confusion() * 100.0
    {
        return Ok(None);
    }
    let bounds = ogeom_algo::shape_bounds(model, profile, tol)?;
    let (Some(low), Some(high)) = (bounds.low(), bounds.high()) else {
        return Ok(None);
    };
    let width = low.distance(high);
    let unit = |v: Vector| v / v.magnitude();
    // Each corner's mitre: its point and the normal of the plane halving it,
    // pointing on along the next leg.
    let mut mitres: Vec<(Point, Vector)> = Vec::with_capacity(legs.len().saturating_sub(1));
    for pair in legs.windows(2) {
        let (d0, d1) = (unit(pair[0].1), unit(pair[1].1));
        // A corner turning past a right angle and a half reaches too far
        // along its legs for the mitre to stay on them.
        if d0.dot(d1) < -0.7 {
            return Ok(None);
        }
        mitres.push((pair[1].0, unit(d0 + d1)));
    }
    // How far a section reaches along its leg from a corner's mitre.
    let reach = width * 2.0;
    for (i, (_, v)) in legs.iter().enumerate() {
        let corners = f64::from(u8::from(i > 0) + u8::from(i + 1 < legs.len()));
        if v.magnitude() <= reach * corners * 0.5 {
            return Ok(None);
        }
    }

    let mut carried = Transform::IDENTITY;
    let mut joined: Vec<Shape> = Vec::new();
    for (i, &(_, v)) in legs.iter().enumerate() {
        let d = unit(v);
        let before = if i > 0 { reach } else { 0.0 };
        let after = if i + 1 < legs.len() { reach } else { 0.0 };
        let section = realized_profile_wound(
            model,
            profile,
            &(Transform::translation(-d * before) * carried),
            Some(d),
            tol,
        )?;
        let mut piece =
            ogeom_algo::make_prism(model, &section, v + d * (before + after), tol)?.shape;
        // Trimmed back to the mitres at either end.
        for (at, normal, keep_ahead) in [
            (i > 0).then(|| (mitres[i - 1].0, mitres[i - 1].1, true)),
            (i + 1 < legs.len()).then(|| (mitres[i].0, mitres[i].1, false)),
        ]
        .into_iter()
        .flatten()
        {
            let plane = Plane::through(at, Direction::new(normal, tol)?);
            let face = ogeom_algo::make_natural_face(
                model,
                SurfaceGeometry::Plane(ogeom_geom::PlaneSurface::over(
                    plane,
                    (-reach * 4.0 - width, reach * 4.0 + width),
                    (-reach * 4.0 - width, reach * 4.0 + width),
                )?),
            )?
            .shape;
            let side = if keep_ahead { normal } else { -normal };
            let half = ogeom_algo::make_half_space(model, &face, at + side * width, tol)?.shape;
            piece = ogeom_bool::common(model, &piece, &half, tol)?.shape;
        }
        // On to the next leg: along this one, then turned at the corner by
        // the least rotation between the headings.
        carried = Transform::translation(v) * carried;
        if let Some(&(corner, next)) = legs.get(i + 1) {
            let n = unit(next);
            let axis = d.cross(n);
            if axis.magnitude() > tol.angular() {
                let angle = d.dot(n).clamp(-1.0, 1.0).acos();
                let rotation = Transform::rotation(
                    ogeom_math::Axis {
                        location: corner,
                        direction: Direction::new(axis, tol)?,
                    },
                    angle,
                );
                carried = rotation * carried;
            }
        }
        joined.push(piece);
    }
    let result = fuse_in_order(model, joined, tol)?;
    Ok(result.map(|shape| {
        let mut history = History::new();
        history.generate(spine, shape.clone());
        history.generate(profile, shape.clone());
        Built::new(shape, history)
    }))
}

/// Sweep a planar profile (a wire, or a face whose holes ride along)
/// down an arbitrary spine, one skinned wall per profile loop.
///
/// The spine may be a single edge or a wire of edges of any curve the
/// vocabulary evaluates: lines, arcs, splines, helices. Frames along it are
/// rotation-minimizing by default (the double-reflection construction), so
/// the profile neither twists nor kinks where the spine bends; `frenet`
/// asks for the Frenet frame instead, which turns with the spine's own
/// curvature, the law a thread wants. Stations are placed by each edge's
/// own turning, the skin holds every transported section to `tolerance`,
/// and the caps sit perpendicular to the spine's ends, holes and all.
///
/// Each spine edge skins its own run of wall, and neighbouring runs share
/// the section where their edges meet, so a join where the curvature steps
/// (an arc running on into its tangent line) is followed exactly rather
/// than smoothed by one fit across it. A run whose sections all lie in one
/// plane is that plane.
///
/// A sharp corner is mitred. Between straight legs the mitre is a plane and
/// each wall is sheared onto it; where a leg is curved the two legs' walls
/// end on the crossing of their generators (each profile point's own path
/// down either leg, run straight on past the corner), which is exact where
/// the corner turns in the leg's plane. Where it turns a curved leg out of
/// its plane the generators miss, and the spine is swept in pieces instead:
/// each side runs on straight past the corner, is trimmed by the mitre
/// plane, and the pieces are fused, the difference between their sections
/// standing as a face of the mitre plane.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the
/// profile is not planar, leans along the spine, or does not sit at the
/// spine's start; if `frenet` is asked of a spine that never bends, or of a
/// cornered one; if a wire (not a face) is swept round a corner that turns
/// a curved leg out of its plane, which only solid pieces can mitre; if the
/// spine all but doubles back at a corner; or if a leg is shorter than its
/// corner's reach.
/// [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) if the skin
/// cannot reach the tolerance.
pub fn make_pipe_shell(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    frenet: bool,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    let law = if frenet {
        PipeLaw::Frenet
    } else {
        PipeLaw::RotationMinimizing
    };
    make_pipe_shell_with(
        model,
        profile,
        spine,
        &law,
        PipeCorners::Mitre,
        tolerance,
        tol,
    )
}

/// How a pipe's section turns about its spine.
#[derive(Debug, Clone, Copy)]
pub enum PipeLaw<'a> {
    /// Double-reflection rotation-minimizing frames: the section neither
    /// twists nor kinks where the spine bends.
    RotationMinimizing,
    /// The Frenet frame, turning with the spine's own curvature.
    Frenet,
    /// The section's normal axis points at `guide`, a curve running beside
    /// the spine: at each station, where the guide crosses the plane square
    /// to the spine there.
    Auxiliary {
        /// An edge or wire beside the spine.
        guide: &'a Shape,
    },
    /// The section keeps this direction as its binormal, square to the
    /// spine's tangent.
    Binormal(Direction),
    /// The section keeps the frame it has at the spine's start all the
    /// way: carried by translation along the spine, never turned.
    Fixed,
}

/// How a pipe turns a sharp corner of its spine.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
pub enum PipeCorners {
    /// The walls sheared onto the plane halving the corner.
    #[default]
    Mitre,
    /// Each leg runs on straight past the corner until its section clears
    /// the other leg's far side, and the legs are fused: the outside of the
    /// corner square.
    Extended,
    /// Each leg ends square at the corner and the section turns about the
    /// corner through its angle, joining the two ends: the outside of the
    /// corner rounded.
    Round,
}

/// [`make_pipe_shell`] with a frame law: the section turns about the spine
/// the way `law` says.
///
/// # Errors
///
/// As [`make_pipe_shell_with`].
pub fn make_pipe_shell_law(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    law: PipeLaw<'_>,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    make_pipe_shell_with(
        model,
        profile,
        spine,
        &law,
        PipeCorners::Mitre,
        tolerance,
        tol,
    )
}

/// [`make_pipe_shell`] with a frame law and a way of turning corners.
///
/// The auxiliary and binormal laws sweep an open spine with no sharp
/// corner. Extended and round corners are built for a face swept down an
/// open spine of straight legs with the rotation-minimizing frame, each
/// leg a prism of the section and each round corner the section revolved
/// about the corner; on a spine with no sharp corner every way of turning
/// corners is the same solid.
///
/// # Errors
///
/// As [`make_pipe_shell`], and
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) where
/// the law or the corners are asked of a spine they do not serve, the guide
/// does not cross a station's plane, or the spine's tangent runs along the
/// binormal or toward the guide.
#[allow(clippy::too_many_arguments)]
pub fn make_pipe_shell_with(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    law: &PipeLaw<'_>,
    corners: PipeCorners,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    // A section carried by translation turns no corner: every way of
    // turning one is the same solid.
    if matches!(law, PipeLaw::Fixed) {
        return crate::fixed::fixed_pipe(model, profile, spine, tol);
    }
    if corners != PipeCorners::Mitre
        && let Some(legs) = straight_legs(model, spine, tol)?
        && legs.len() > 1
    {
        if !matches!(law, PipeLaw::RotationMinimizing) {
            ogeom_bail!(
                Construction,
                "extended and round corners turn the section with the \
                 rotation-minimizing frame"
            );
        }
        if model.kind_of(profile)? != ShapeType::Face {
            ogeom_bail!(Construction, "extended and round corners sweep a face");
        }
        return cornered_lines(model, profile, spine, &legs, corners, tol);
    }
    if corners != PipeCorners::Mitre && has_sharp_corner(model, spine, tol)? {
        ogeom_bail!(
            Construction,
            "extended and round corners are built along a spine of straight legs"
        );
    }
    pipe_shell_law(model, profile, spine, law, tolerance, tol)
}

/// The legs of a spine made only of straight edges, in order, each as its
/// start and its displacement; `None` where any edge is curved.
fn straight_legs(
    model: &Model,
    spine: &Shape,
    tol: Tolerances,
) -> OgeomResult<Option<Vec<(Point, Vector)>>> {
    use ogeom_geom::Curve;
    let edges: Vec<Shape> = match model.kind_of(spine)? {
        ShapeType::Edge => vec![spine.clone()],
        ShapeType::Wire => model.ordered_children_of(spine)?,
        _ => return Ok(None),
    };
    let mut legs = Vec::with_capacity(edges.len());
    for edge in &edges {
        let (curve, range) = spine_curve_of(model, edge)?;
        let basis = match &curve {
            Curve::Trimmed(t) => t.basis().clone(),
            other => other.clone(),
        };
        if !matches!(basis, Curve::Line(_)) {
            return Ok(None);
        }
        let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
        let (t0, t1) = if reversed { (range.1, range.0) } else { range };
        let (a, b) = (curve.point_at(t0, tol)?, curve.point_at(t1, tol)?);
        legs.push((a, b - a));
    }
    Ok(Some(legs))
}

/// Whether the spine turns sharply anywhere between its edges.
pub(crate) fn has_sharp_corner(model: &Model, spine: &Shape, tol: Tolerances) -> OgeomResult<bool> {
    let stations = shell_stations(model, spine, tol)?;
    Ok((0..stations.len() - 1).any(|i| {
        stations[i].at.distance(stations[i + 1].at) <= tol.confusion()
            && (stations[i]
                .tangent
                .cross(stations[i + 1].tangent)
                .magnitude()
                > tol.angular()
                || stations[i].tangent.dot(stations[i + 1].tangent) < 0.0)
    }))
}

/// A face down an open spine of straight legs with extended or round
/// corners: each leg the prism of its section, run on past its corners
/// for extended ones, and each round corner the section revolved about the
/// corner through its turn; the pieces fused.
fn cornered_lines(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    legs: &[(Point, Vector)],
    corners: PipeCorners,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if let (Some(first), Some(last)) = (legs.first(), legs.last())
        && first.0.distance(last.0 + last.1) <= tol.confusion() * 100.0
    {
        ogeom_bail!(
            Construction,
            "extended and round corners are built along an open spine"
        );
    }
    let unit = |v: Vector| v / v.magnitude();
    let Some(outer) = explore(model, profile, Filter::OfType(ShapeType::Wire))?
        .into_iter()
        .next()
    else {
        ogeom_bail!(Construction, "the profile has no loop");
    };
    let samples = sample_wire(model, &outer, 256, tol)?;
    // The motion carrying the profile to each leg's start: along the legs
    // before it, turned at each corner by the least rotation between the
    // headings.
    let mut carried: Vec<Transform> = Vec::with_capacity(legs.len());
    let mut turns: Vec<Option<(ogeom_math::Axis, f64)>> = Vec::with_capacity(legs.len());
    let mut motion = Transform::IDENTITY;
    for (i, &(_, v)) in legs.iter().enumerate() {
        carried.push(motion);
        motion = Transform::translation(v) * motion;
        let turn = match legs.get(i + 1) {
            Some(&(corner, next)) => {
                let (d, n) = (unit(v), unit(next));
                if d.dot(n) < -0.7 {
                    ogeom_bail!(
                        Construction,
                        "the spine all but doubles back at a corner; no corner \
                         of that kind turns it"
                    );
                }
                let axis = d.cross(n);
                if axis.magnitude() > tol.angular() {
                    let about = ogeom_math::Axis {
                        location: corner,
                        direction: Direction::new(axis, tol)?,
                    };
                    let angle = d.dot(n).clamp(-1.0, 1.0).acos();
                    motion = Transform::rotation(about, angle) * motion;
                    Some((about, angle))
                } else {
                    None
                }
            }
            None => None,
        };
        turns.push(turn);
    }
    // How far each leg runs on past a corner for an extended corner: until
    // its section passes the far side of the other leg's section there.
    let reach_past = |section: &Transform, corner: Point, along: Vector, turn: f64| {
        let far = samples
            .iter()
            .map(|p| (section.apply(*p) - corner).dot(along))
            .fold(0.0_f64, f64::max);
        far / turn.sin().max(1e-3)
    };
    let mut joined: Vec<Shape> = Vec::new();
    for (i, &(_, v)) in legs.iter().enumerate() {
        let d = unit(v);
        let (mut before, mut after) = (0.0, 0.0);
        if corners == PipeCorners::Extended {
            if i > 0
                && let Some((_, angle)) = turns[i - 1]
            {
                // The previous leg's section at the corner, reaching back.
                let at_corner = Transform::translation(legs[i - 1].1) * carried[i - 1];
                before = reach_past(&at_corner, legs[i].0, -d, angle);
            }
            if let Some((_, angle)) = turns[i] {
                // The next leg's section at the corner, reaching on.
                after = reach_past(&carried[i + 1], legs[i + 1].0, d, angle);
            }
        }
        let section = realized_profile_wound(
            model,
            profile,
            &(Transform::translation(-d * before) * carried[i]),
            Some(d),
            tol,
        )?;
        let mut piece =
            ogeom_algo::make_prism(model, &section, v + d * (before + after), tol)?.shape;
        if corners == PipeCorners::Round
            && let Some((axis, angle)) = turns[i]
        {
            let end = realized_profile_wound(
                model,
                profile,
                &(Transform::translation(v) * carried[i]),
                Some(d),
                tol,
            )?;
            // The half of the section inside the turn sweeps within the two
            // legs; the half outside it turns the corner's round. It
            // touches the axis along its cut, which the revolution takes.
            let next = unit(legs[i + 1].1);
            let inward = next - d * next.dot(d);
            let outer = outside_half(model, &end, axis.location, d, -inward, tol)?;
            let bend = ogeom_algo::make_revolution(model, &outer, axis, angle, tol)?.shape;
            piece = ogeom_bool::fuse(model, &piece, &bend, tol)?.shape;
        }
        joined.push(piece);
    }
    let result = fuse_in_order(model, joined, tol)?;
    let Some(shape) = result else {
        ogeom_bail!(Construction, "the spine has no leg");
    };
    let mut history = History::new();
    history.generate(spine, shape.clone());
    history.generate(profile, shape.clone());
    Ok(Built::new(shape, history))
}

/// The part of a planar section lying on the `outward` side of the line
/// through `at` square to it in the section's plane (whose normal is
/// `normal`), as a face on that plane: the section's slab kept on that
/// side by a half space, and the slab's face on the section's plane.
fn outside_half(
    model: &mut Model,
    section: &Shape,
    at: Point,
    normal: Vector,
    outward: Vector,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    let width = {
        let bounds = ogeom_algo::shape_bounds(model, section, tol)?;
        match (bounds.low(), bounds.high()) {
            (Some(a), Some(b)) => a.distance(b),
            _ => ogeom_bail!(Construction, "the section has no extent"),
        }
    };
    let n = normal / normal.magnitude();
    let slab = ogeom_algo::make_prism(model, section, n * width, tol)?.shape;
    let cut = Plane::through(at, Direction::new(outward, tol)?);
    let face = ogeom_algo::make_natural_face(
        model,
        SurfaceGeometry::Plane(ogeom_geom::PlaneSurface::over(
            cut,
            (-width * 4.0, width * 4.0),
            (-width * 4.0, width * 4.0),
        )?),
    )?
    .shape;
    let side = outward / outward.magnitude() * width;
    let half = ogeom_algo::make_half_space(model, &face, at + side, tol)?.shape;
    let kept = ogeom_bool::common(model, &slab, &half, tol)?.shape;
    for face in explore(model, &kept, Filter::OfType(ShapeType::Face))? {
        let Ok((point, facing)) = ogeom_algo::face_normal(model, &face, tol) else {
            continue;
        };
        if (point - at).dot(n).abs() <= tol.confusion() * 100.0
            && facing.cross(n).magnitude() <= tol.angular() * 10.0
        {
            // Facing along the section's own normal, as the section does.
            return Ok(if facing.dot(n) > 0.0 {
                face
            } else {
                face.reversed()
            });
        }
    }
    ogeom_bail!(
        Construction,
        "the section has no part outside the corner's axis to turn"
    )
}

/// Stations added between neighbours on one edge wherever a law's frame
/// turns by more than a few degrees between them: a straight spine has
/// stations only at its ends, and a section twisting down it would be
/// skinned straight across.
fn densified(
    model: &Model,
    spine: &Shape,
    stations: Vec<SpineStation>,
    law: &PipeLaw<'_>,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Vec<SpineStation>> {
    const MOST: f64 = 0.05;
    let edges: Vec<Shape> = match model.kind_of(spine)? {
        ShapeType::Edge => vec![spine.clone()],
        ShapeType::Wire => model.ordered_children_of(spine)?,
        _ => return Ok(stations),
    };
    let mut curves = Vec::with_capacity(edges.len());
    for edge in &edges {
        curves.push(spine_curve_of(model, edge)?.0);
    }
    let mut stations = stations;
    for _ in 0..12 {
        let normals = law_normals(model, &stations, law, tolerance, tol)?;
        let mut out: Vec<SpineStation> = Vec::with_capacity(stations.len() * 2);
        let mut added = false;
        for i in 0..stations.len() {
            out.push(stations[i]);
            let Some(next) = stations.get(i + 1) else {
                continue;
            };
            let here = stations[i];
            if here.edge != next.edge || here.at.distance(next.at) <= tol.confusion() {
                continue;
            }
            let turn = normals[i].dot(normals[i + 1]).clamp(-1.0, 1.0).acos();
            if turn <= MOST {
                continue;
            }
            let curve = &curves[here.edge];
            let t = f64::midpoint(here.t, next.t);
            let d = curve.d1_at(t, tol)?;
            let sense = if curve.d1_at(here.t, tol)?.dot(here.tangent) >= 0.0 {
                1.0
            } else {
                -1.0
            };
            out.push(SpineStation {
                at: curve.point_at(t, tol)?,
                tangent: d * sense / d.magnitude(),
                edge: here.edge,
                t,
            });
            added = true;
        }
        stations = out;
        if !added {
            break;
        }
        // A law whose frame jumps (a guide's nearest point switching
        // branch) never settles: each round doubles the stations at the
        // jump. Past this many the skin is refused, not grown.
        if stations.len() > MOST_STATIONS {
            ogeom_bail!(
                NotDone,
                "the frame law turns faster than {MOST_STATIONS} sections can follow; \
                 it jumps somewhere along the spine"
            );
        }
    }
    Ok(stations)
}

/// The most sections a law pipe's skin is built through.
const MOST_STATIONS: usize = 2048;

/// The frame a law gives at a station: the spine's tangent as its `z`, the
/// law's normal as its `x`.
pub(crate) fn station_frame(
    station: &SpineStation,
    normal: Vector,
    tol: Tolerances,
) -> OgeomResult<Frame> {
    Frame::new(
        station.at,
        Direction::new(station.tangent, tol)?,
        Direction::new(normal, tol)?,
        tol,
    )
}

/// A copy of a closed wire moved by `motion`, its curves restated where
/// they land rather than placed: every edge in the wire's order and sense,
/// on vertices shared end to end.
fn moved_ring(
    model: &mut Model,
    ring: &Shape,
    motion: &Transform,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    let edges = model.ordered_children_of(ring)?;
    let mut vertices: Vec<(Point, Shape)> = Vec::new();
    let mut vertex = |model: &mut Model, at: Point| -> Shape {
        if let Some((_, v)) = vertices
            .iter()
            .find(|(p, _)| p.distance(at) <= tol.confusion() * 10.0)
        {
            return v.clone();
        }
        let v = ogeom_algo::make_vertex(model, at).shape;
        vertices.push((at, v.clone()));
        v
    };
    let mut moved = Vec::with_capacity(edges.len());
    for edge in &edges {
        let (curve, range) = spine_curve_of(model, edge)?;
        let placement = edge.transform(model.datums())?;
        let curve = curve.transformed(&(*motion * placement), tol)?;
        let (a, b) = (curve.point_at(range.0, tol)?, curve.point_at(range.1, tol)?);
        let (va, vb) = (vertex(model, a), vertex(model, b));
        let built = ogeom_algo::make_edge_between(model, curve, range, &va, &vb, tol)?.shape;
        moved.push(if edge.orientation() == ogeom_topo::Orientation::Reversed {
            built.reversed()
        } else {
            built
        });
    }
    Ok(ogeom_algo::make_wire(model, &moved, tol)?.shape)
}

/// A profile swept under a law that turns it about the spine, as the loft
/// through its outer ring placed at every station by the law's frame: the
/// section keeps the place it has in the frame at the spine's start.
#[allow(clippy::too_many_arguments)]
fn law_loft(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    stations: &[SpineStation],
    law: &PipeLaw<'_>,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    let rings: Vec<Shape> = match model.kind_of(profile)? {
        ShapeType::Face => model.ordered_children_of(profile)?,
        ShapeType::Wire => vec![profile.clone()],
        _ => ogeom_bail!(Construction, "a pipe sweeps a planar face or wire"),
    };
    if rings.len() != 1 {
        ogeom_bail!(
            Construction,
            "an auxiliary or binormal law sweeps a profile with no hole"
        );
    }
    let normals = law_normals(model, stations, law, tolerance, tol)?;
    let start = station_frame(&stations[0], normals[0], tol)?;
    let mut sections: Vec<Shape> = Vec::with_capacity(stations.len());
    let mut last: Option<Point> = None;
    for (station, normal) in stations.iter().zip(&normals) {
        if last.is_some_and(|p| p.distance(station.at) <= tol.confusion()) {
            continue;
        }
        last = Some(station.at);
        let frame = station_frame(station, *normal, tol)?;
        let motion = Transform::from_frame(&frame) * Transform::to_frame(&start);
        sections.push(moved_ring(model, &rings[0], &motion, tol)?);
    }
    let mut built = make_loft_skinned(model, &sections, tolerance, tol)?;
    built.history.generate(spine, built.shape.clone());
    built.history.generate(profile, built.shape.clone());
    Ok(built)
}

/// Sweep several planar sections down one spine, the section changing
/// shape along the path: a multisection pipe.
///
/// Each section stands where the spine crosses its plane, and is read in
/// the spine's moving frame there (rotation-minimizing, or Frenet where
/// `frenet` asks). Between two sections the section is blended in that
/// frame by the length run along the spine, matched point to point from
/// each section's own start, so the result follows the spine rather than
/// the chord between the sections. The sections are closed wires, or
/// faces without holes; the ends are capped. The first or last section may
/// be a vertex where the spine starts or ends: the pipe closes to that
/// point, uncapped there.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if
/// there are fewer than two sections, a section is not planar or has a
/// hole, the spine does not cross a section's plane, or the sections stand
/// out of order along the spine;
/// [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) if the skin
/// cannot reach `tolerance`.
pub fn make_pipe_sections(
    model: &mut Model,
    sections: &[Shape],
    spine: &Shape,
    frenet: bool,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    const AROUND: usize = 64;
    if sections.len() < 2 {
        ogeom_bail!(
            Construction,
            "a multisection pipe needs at least two sections"
        );
    }
    let law = if frenet {
        PipeLaw::Frenet
    } else {
        PipeLaw::RotationMinimizing
    };
    let stations = shell_stations(model, spine, tol)?;
    let stations = evenly(model, spine, stations, 24, false, tol)?;
    let normals = law_normals(model, &stations, &law, tolerance, tol)?;
    // Length run along the spine at each station.
    let mut run = vec![0.0_f64];
    for pair in stations.windows(2) {
        let held = run[run.len() - 1];
        run.push(held + pair[0].at.distance(pair[1].at));
    }
    // Each section: where along the spine it stands, and its samples in
    // the frame there.
    let mut placed: Vec<(f64, Vec<Point>)> = Vec::with_capacity(sections.len());
    // A point section: where it stands, and whether at the spine's end.
    let mut apex: Option<(Point, bool)> = None;
    // The last section's ring read, for the straight pipe to a point.
    let mut last_ring: Option<Shape> = None;
    for (index, section) in sections.iter().enumerate() {
        if model.kind_of(section)? == ShapeType::Vertex {
            let Some(data) = model.node(section).and_then(|n| n.data().as_vertex()) else {
                ogeom_bail!(Dangling, "vertex is not in this model");
            };
            let point = section.transform(model.datums())?.apply(data.point);
            let reach = tolerance.max(tol.confusion() * 1e3);
            let (Some(head), Some(tail)) = (stations.first(), stations.last()) else {
                ogeom_bail!(Construction, "the spine has no stations");
            };
            let at_end = if index + 1 == sections.len() && point.distance(tail.at) <= reach {
                true
            } else if index == 0 && point.distance(head.at) <= reach {
                false
            } else {
                ogeom_bail!(
                    Construction,
                    "a point section stands first or last, where the spine \
                     starts or ends"
                );
            };
            apex = Some((point, at_end));
            // Every sample at the frame's origin: the section shrinks onto
            // the spine there.
            let along = if at_end { run[run.len() - 1] } else { 0.0 };
            placed.push((along, vec![Point::ORIGIN; AROUND]));
            continue;
        }
        let ring = match model.kind_of(section)? {
            ShapeType::Face => {
                let rings = model.ordered_children_of(section)?;
                if rings.len() != 1 {
                    ogeom_bail!(Construction, "a multisection pipe's section has no hole");
                }
                rings[0].clone()
            }
            ShapeType::Wire => section.clone(),
            _ => ogeom_bail!(
                Construction,
                "a section is a planar face or wire, or a point at an end"
            ),
        };
        let Some(plane) = ogeom_algo::find_plane(model, &ring, tol)? else {
            ogeom_bail!(Construction, "a section is not planar");
        };
        last_ring = Some(ring.clone());
        let side = |p: Point| (p - plane.origin()).dot(plane.normal().vector());
        let mut found: Option<(f64, usize, f64)> = None;
        for i in 0..stations.len() {
            let here = side(stations[i].at);
            let crossing = if here.abs() <= tol.confusion() {
                Some((run[i], i, 0.0))
            } else if let Some(next) = stations.get(i + 1) {
                let there = side(next.at);
                (here.signum() != there.signum() && there.abs() > tol.confusion()).then(|| {
                    let f = here / (here - there);
                    (run[i] + (run[i + 1] - run[i]) * f, i, f)
                })
            } else {
                None
            };
            if crossing.is_some() {
                found = crossing;
                break;
            }
        }
        let Some((along, i, f)) = found else {
            ogeom_bail!(Construction, "the spine does not cross a section's plane");
        };
        // The frame there, between the two stations it falls between.
        let j = (i + 1).min(stations.len() - 1);
        let at = stations[i].at + (stations[j].at - stations[i].at) * f;
        let tangent = stations[i].tangent * (1.0 - f) + stations[j].tangent * f;
        let normal = normals[i] * (1.0 - f) + normals[j] * f;
        let normal = normal - tangent * normal.dot(tangent) / tangent.dot(tangent);
        let frame = Frame::new(
            at,
            Direction::new(tangent, tol)?,
            Direction::new(normal, tol)?,
            tol,
        )?;
        let into = Transform::to_frame(&frame);
        let dense: Vec<Point> = dense_wire(model, &ring, tol)?
            .iter()
            .map(|p| into.apply(*p))
            .collect();
        let samples = match placed.last() {
            Some((_, previous)) => matched_samples(&dense, AROUND, previous),
            None => resampled(&dense, AROUND),
        };
        placed.push((along, samples));
    }
    for pair in placed.windows(2) {
        if pair[1].0 <= pair[0].0 + tol.confusion() {
            ogeom_bail!(
                Construction,
                "the sections stand out of order along the spine"
            );
        }
    }
    // One section closing to a point down a straight spine: every sample
    // runs straight to the point, which is the exact pyramid or cone over
    // the section.
    if let (Some((point, _)), Some(ring), 2) = (apex, &last_ring, sections.len())
        && straight_legs(model, spine, tol)?.is_some_and(|legs| legs.len() == 1)
    {
        let tip = model.add_vertex(VertexData::new(point));
        let mut built = make_loft(model, ring, &tip, tol)?;
        built.history.generate(spine, built.shape.clone());
        for section in sections {
            built.history.generate(section, built.shape.clone());
        }
        return Ok(built);
    }
    // A ring at every station from the first section's to the last's, the
    // two sections either side blended by the length between them.
    let (first, last) = (placed[0].0, placed[placed.len() - 1].0);
    let mut rings: Vec<Shape> = Vec::new();
    let mut ring_at = |model: &mut Model, along: f64, frame: &Frame| -> OgeomResult<()> {
        let k = placed
            .windows(2)
            .position(|w| along <= w[1].0 + tol.confusion())
            .unwrap_or(placed.len() - 2);
        let (a, b) = (&placed[k], &placed[k + 1]);
        let f = ((along - a.0) / (b.0 - a.0)).clamp(0.0, 1.0);
        let out = Transform::from_frame(frame);
        let mut points: Vec<Point> =
            a.1.iter()
                .zip(&b.1)
                .map(|(p, q)| out.apply(*p + (*q - *p) * f))
                .collect();
        points.push(points[0]);
        let fitted = ogeom_geom::fit::fit_points_closed(&points, 3, tolerance * 0.1, tol)?;
        let curve: ogeom_geom::Curve = fitted.curve.into();
        let range = curve.domain();
        let edge = ogeom_algo::make_edge(model, curve, range, tol)?.shape;
        rings.push(ogeom_algo::make_wire(model, &[edge], tol)?.shape);
        Ok(())
    };
    let frame_at = |i: usize| station_frame(&stations[i], normals[i], tol);
    // The first section's own station, the stations strictly between, and
    // the last section's.
    let at_along = |along: f64| -> OgeomResult<Frame> {
        let i = run
            .windows(2)
            .position(|w| along <= w[1] + tol.confusion())
            .unwrap_or(run.len() - 2);
        let f = ((along - run[i]) / (run[i + 1] - run[i]).max(f64::MIN_POSITIVE)).clamp(0.0, 1.0);
        let at = stations[i].at + (stations[i + 1].at - stations[i].at) * f;
        let tangent = stations[i].tangent * (1.0 - f) + stations[i + 1].tangent * f;
        let normal = normals[i] * (1.0 - f) + normals[i + 1] * f;
        let normal = normal - tangent * normal.dot(tangent) / tangent.dot(tangent);
        Frame::new(
            at,
            Direction::new(tangent, tol)?,
            Direction::new(normal, tol)?,
            tol,
        )
    };
    // The end a point section stands at has no ring: the skin closes on the
    // point itself, which the loft takes last.
    if apex.is_none_or(|(_, at_end)| at_end) {
        ring_at(model, first, &at_along(first)?)?;
    }
    for (i, &along) in run.iter().enumerate() {
        if along > first + tol.confusion() && along < last - tol.confusion() {
            ring_at(model, along, &frame_at(i)?)?;
        }
    }
    if apex.is_none_or(|(_, at_end)| !at_end) {
        ring_at(model, last, &at_along(last)?)?;
    }
    if let Some((point, at_end)) = apex {
        if !at_end {
            rings.reverse();
        }
        rings.push(model.add_vertex(VertexData::new(point)));
    }
    let mut built = make_loft_skinned(model, &rings, tolerance, tol)?;
    built.history.generate(spine, built.shape.clone());
    for section in sections {
        built.history.generate(section, built.shape.clone());
    }
    Ok(built)
}

/// Stations added so every edge of the spine holds at least `count`
/// stations evenly along its parameter, besides the ones it has.
fn evenly(
    model: &Model,
    spine: &Shape,
    stations: Vec<SpineStation>,
    count: u32,
    keep: bool,
    tol: Tolerances,
) -> OgeomResult<Vec<SpineStation>> {
    evenly_by(model, spine, stations, &|_| count, keep, tol)
}

/// [`evenly`], with each edge's own count: `count(e)` for edge `e`.
fn evenly_by(
    model: &Model,
    spine: &Shape,
    stations: Vec<SpineStation>,
    count: &dyn Fn(usize) -> u32,
    keep: bool,
    tol: Tolerances,
) -> OgeomResult<Vec<SpineStation>> {
    let edges: Vec<Shape> = match model.kind_of(spine)? {
        ShapeType::Edge => vec![spine.clone()],
        ShapeType::Wire => model.ordered_children_of(spine)?,
        _ => return Ok(stations),
    };
    let mut out: Vec<SpineStation> = Vec::with_capacity(stations.len() + edges.len() * 32);
    for (e, edge) in edges.iter().enumerate() {
        let on: Vec<SpineStation> = stations.iter().copied().filter(|s| s.edge == e).collect();
        let (Some(first), Some(last)) = (on.first().copied(), on.last().copied()) else {
            continue;
        };
        let (curve, _) = spine_curve_of(model, edge)?;
        let sense = if curve.d1_at(first.t, tol)?.dot(first.tangent) >= 0.0 {
            1.0
        } else {
            -1.0
        };
        let mut ts: Vec<f64> = if keep {
            on.iter().map(|s| s.t).collect()
        } else {
            vec![first.t, last.t]
        };
        let count = count(e);
        for k in 1..count {
            ts.push(first.t + (last.t - first.t) * f64::from(k) / f64::from(count));
        }
        ts.sort_by(|a, b| {
            if first.t <= last.t {
                a.total_cmp(b)
            } else {
                b.total_cmp(a)
            }
        });
        ts.dedup_by(|a, b| (*a - *b).abs() <= tol.parametric());
        for t in ts {
            let d = curve.d1_at(t, tol)?;
            out.push(SpineStation {
                at: curve.point_at(t, tol)?,
                tangent: d * sense / d.magnitude(),
                edge: e,
                t,
            });
        }
    }
    Ok(out)
}

/// The frame normals a law gives at each station.
///
/// An auxiliary guide ending within `reach` of a station's plane (a
/// sketch's end in single precision) is taken to cross it, carried on
/// along its end tangent.
pub(crate) fn law_normals(
    model: &Model,
    stations: &[SpineStation],
    law: &PipeLaw<'_>,
    reach: f64,
    tol: Tolerances,
) -> OgeomResult<Vec<Vector>> {
    match law {
        // Carried by translation, the section has no frame turning with the
        // spine; the fixed pipe is built without stations.
        PipeLaw::Fixed => ogeom_bail!(
            Construction,
            "a fixed section keeps its own frame; it has no frame normals"
        ),
        PipeLaw::RotationMinimizing => Ok(rmf_normals(stations)),
        PipeLaw::Frenet => frenet_normals(stations, tol),
        PipeLaw::Binormal(b) => stations
            .iter()
            .map(|s| {
                let t = s.tangent;
                let b = b.vector() - t * b.vector().dot(t);
                if b.magnitude() <= tol.angular() {
                    ogeom_bail!(
                        Construction,
                        "the spine runs along the binormal at {:?}; no frame keeps it",
                        s.at
                    );
                }
                let n = (b / b.magnitude()).cross(t);
                Ok(n / n.magnitude())
            })
            .collect(),
        PipeLaw::Auxiliary { guide } => {
            let edges: Vec<Shape> = match model.kind_of(guide)? {
                ShapeType::Edge => vec![(*guide).clone()],
                ShapeType::Wire => model.ordered_children_of(guide)?,
                _ => ogeom_bail!(Construction, "an auxiliary spine is an edge or a wire"),
            };
            let mut curves = Vec::with_capacity(edges.len());
            for edge in &edges {
                curves.push(spine_curve_of(model, edge)?);
            }
            let mut out = Vec::with_capacity(stations.len());
            let mut last: Option<Point> = None;
            for s in stations {
                let (p, t) = (s.at, s.tangent);
                // Where the guide crosses this station's plane: a sign
                // change of the height along each guide edge, refined; the
                // crossing nearest the last one found.
                let mut best: Option<Point> = None;
                for (curve, range) in &curves {
                    const STEPS: u32 = 256;
                    let height = |u: f64| -> OgeomResult<(f64, Point)> {
                        let q = curve.point_at(u, tol)?;
                        Ok(((q - p).dot(t), q))
                    };
                    let at =
                        |k: u32| range.0 + (range.1 - range.0) * f64::from(k) / f64::from(STEPS);
                    let mut prev = height(at(0))?;
                    for k in 1..=STEPS {
                        let here = height(at(k))?;
                        if prev.0 == 0.0 || prev.0.signum() != here.0.signum() {
                            let (mut lo, mut hi) = (at(k - 1), at(k));
                            let mut f_lo = prev.0;
                            for _ in 0..60 {
                                let mid = f64::midpoint(lo, hi);
                                let (f_mid, _) = height(mid)?;
                                if f_mid.signum() == f_lo.signum() {
                                    lo = mid;
                                    f_lo = f_mid;
                                } else {
                                    hi = mid;
                                }
                            }
                            let q = height(f64::midpoint(lo, hi))?.1;
                            let near = last.unwrap_or(p);
                            if best.is_none_or(|b| q.distance(near) < b.distance(near)) {
                                best = Some(q);
                            }
                        }
                        prev = here;
                    }
                }
                if best.is_none() {
                    let near = last.unwrap_or(p);
                    for (curve, range) in &curves {
                        for u in [range.0, range.1] {
                            let q = curve.point_at(u, tol)?;
                            let h = (q - p).dot(t);
                            if h.abs() > reach {
                                continue;
                            }
                            let d = curve.d1_at(u, tol)?;
                            let along = d.dot(t);
                            let onto = if along.abs() > tol.angular() * d.magnitude() {
                                q - d * (h / along)
                            } else {
                                q
                            };
                            // A tangent nearly in the plane carries the end
                            // far; the end itself is then as good.
                            let q = if onto.distance(q) <= 2.0 * reach {
                                onto
                            } else {
                                q
                            };
                            if best.is_none_or(|b| q.distance(near) < b.distance(near)) {
                                best = Some(q);
                            }
                        }
                    }
                }
                let Some(q) = best else {
                    ogeom_bail!(
                        Construction,
                        "the auxiliary spine does not cross the plane square to the \
                         spine at {p:?}"
                    );
                };
                last = Some(q);
                let toward = (q - p) - t * (q - p).dot(t);
                if toward.magnitude() <= tol.confusion() {
                    ogeom_bail!(
                        Construction,
                        "the auxiliary spine meets the spine at {p:?}; no direction \
                         points at it"
                    );
                }
                out.push(toward / toward.magnitude());
            }
            Ok(out)
        }
    }
}

/// The pipe shell under a frame law, mitred at corners.
#[allow(clippy::too_many_lines)]
fn pipe_shell_law(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    law: &PipeLaw<'_>,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    const AROUND: usize = 40;
    let frenet = matches!(law, PipeLaw::Frenet);
    let classic = matches!(law, PipeLaw::RotationMinimizing | PipeLaw::Frenet);

    if classic && let Some(exact) = exact_legs(model, profile, spine, frenet, tol)? {
        return Ok(exact);
    }
    let stations = shell_stations(model, spine, tol)?;
    let stations = if classic {
        stations
    } else {
        // As many stations as the law's turning asks, spread evenly: the
        // skin takes its sections at even steps of its own parameter.
        let wanted = densified(model, spine, stations.clone(), law, tolerance, tol)?.len();
        // Shared out along the spine by length, so the spacing is the same
        // on every edge: the skin steps evenly along its whole run.
        let mut lengths: Vec<f64> = Vec::new();
        for pair in stations.windows(2) {
            let e = pair[1].edge;
            if lengths.len() <= e {
                lengths.resize(e + 1, 0.0);
            }
            if pair[0].edge == e {
                lengths[e] += pair[0].at.distance(pair[1].at);
            }
        }
        let total: f64 = lengths.iter().sum();
        #[allow(
            clippy::cast_precision_loss,
            clippy::cast_possible_truncation,
            clippy::cast_sign_loss
        )]
        let share = |e: usize| -> u32 {
            let length = lengths.get(e).copied().unwrap_or(0.0);
            if total <= 0.0 {
                return 2;
            }
            ((wanted as f64 * length / total).ceil() as u32).max(2)
        };
        evenly_by(model, spine, stations, &share, false, tol)?
    };
    // Corners: twin stations standing on one point with different headings.
    let corners: Vec<usize> = (0..stations.len() - 1)
        .filter(|&i| {
            stations[i].at.distance(stations[i + 1].at) <= tol.confusion()
                && (stations[i]
                    .tangent
                    .cross(stations[i + 1].tangent)
                    .magnitude()
                    > tol.angular()
                    || stations[i].tangent.dot(stations[i + 1].tangent) < 0.0)
        })
        .collect();
    // Every twin, corners and smooth junctions alike: each ends one run of
    // skin and starts the next, the two runs sharing the section there. A
    // smooth junction's mitre plane is its own section, so nothing shears.
    let kinks: Vec<usize> = (0..stations.len() - 1)
        .filter(|&i| stations[i].at.distance(stations[i + 1].at) <= tol.confusion())
        .collect();
    let ring = stations[0].at.distance(stations[stations.len() - 1].at) <= tol.confusion() * 10.0;
    if !classic && (ring || !corners.is_empty()) {
        ogeom_bail!(
            Construction,
            "an auxiliary or binormal law sweeps an open spine with no sharp corner"
        );
    }
    if !classic {
        return law_loft(model, profile, spine, &stations, law, tolerance, tol);
    }
    if ring && kinks.is_empty() {
        return closed_pipe_shell(model, profile, spine, stations, frenet, tolerance, tol);
    }
    // The cornered ring closes at its own wrap. Seamed on a corner, the
    // wrap is one more mitre; seamed mid-leg, the wrap's "mitre" plane is
    // the leg's own cross-section (both twin tangents are the leg's), and
    // the shear onto it moves nothing, so the two halves of that leg butt
    // together on the seam's own ring, coplanar walls meeting on it. The
    // solid is exact either way; the mid-leg seam merely leaves its leg in
    // two pieces.
    if frenet && !corners.is_empty() {
        ogeom_bail!(
            Construction,
            "a Frenet frame has no direction at a corner; sweep a cornered \
             spine with the rotation-minimizing frame"
        );
    }
    let normals = law_normals(model, &stations, law, tolerance, tol)?;
    // A ring's frame must come home: carry once more across the wrap
    // corner, read the twist between departure and return, and spread it
    // along the arc: the smooth loop's own reconciliation, ending at a
    // mitre instead of a tangent join.
    let normals = if ring {
        let mut extended = stations.clone();
        extended.push(stations[0]);
        let carried = rmf_normals(&extended);
        let (n0, n_home) = (carried[0], carried[carried.len() - 1]);
        let t0 = stations[0].tangent;
        let twist = (n0.cross(n_home).dot(t0)).atan2(n0.dot(n_home));
        let mut lengths = vec![0.0_f64];
        for pair in extended.windows(2) {
            let held = lengths[lengths.len() - 1];
            lengths.push(held + pair[0].at.distance(pair[1].at));
        }
        let total = lengths[lengths.len() - 1];
        carried
            .iter()
            .take(stations.len())
            .enumerate()
            .map(|(i, n)| {
                let phi = -twist * lengths[i] / total;
                let t = extended[i].tangent;
                let v = *n * phi.cos() + t.cross(*n) * phi.sin();
                let v = v - t * v.dot(t);
                v / v.magnitude()
            })
            .collect()
    } else {
        normals
    };
    // At each corner both twin sections are thrown onto the bisector plane
    // along their own tangents; the mirror symmetry of the rotation-
    // minimizing frame lands them on one ring, the mitre both runs share.
    let mitre: Vec<Option<(Point, Vector)>> = {
        let mut out: Vec<Option<(Point, Vector)>> = vec![None; stations.len()];
        for &k in &kinks {
            let n = stations[k].tangent + stations[k + 1].tangent;
            if n.magnitude() <= tol.angular() {
                ogeom_bail!(
                    Construction,
                    "the spine doubles straight back on itself; no mitre \
                     plane divides that corner"
                );
            }
            out[k] = Some((stations[k].at, n));
            out[k + 1] = Some((stations[k + 1].at, n));
        }
        if ring {
            let wrap = stations.len() - 1;
            let n = stations[wrap].tangent + stations[0].tangent;
            if n.magnitude() <= tol.angular() {
                ogeom_bail!(
                    Construction,
                    "the spine doubles straight back on itself; no mitre \
                     plane divides that corner"
                );
            }
            out[wrap] = Some((stations[wrap].at, n));
            out[0] = Some((stations[0].at, n));
        }
        out
    };
    let runs: Vec<(usize, usize)> = {
        let mut out = Vec::with_capacity(kinks.len() + 1);
        let mut start = 0;
        for &k in &kinks {
            out.push((start, k));
            start = k + 1;
        }
        out.push((start, stations.len() - 1));
        out
    };
    // A mitred end between straight legs is a *shear*: the honest wall is
    // the run's own surface trimmed by the mitre plane, which for a
    // straight leg is exactly the ruled skin between its two end rings.
    // Under a law that turns the section along a straight leg, the leg's
    // wall is no ruled skin between its end rings.
    let straight = |rs: usize, re: usize| -> bool {
        if !classic {
            return false;
        }
        let t0 = stations[rs].tangent;
        (rs..=re).all(|i| stations[i].tangent.cross(t0).magnitude() <= tol.angular())
    };
    // Every corner as the pair of runs it stands between, the wrap
    // included, and whether either leg is curved. A curved leg's trim is
    // not a loft of its rows: the two legs' generators for one profile
    // point are followed (the leg's own curve, run straight on past the
    // corner) to where they meet, and each wall ends on that crossing.
    // Where the corner turns in the plane the crossing is exact; a skew
    // corner's generators miss each other, and that miss is refused.
    struct CornerPair {
        before: (usize, usize),
        after: (usize, usize),
        curved: bool,
    }
    let corner_pairs: Vec<CornerPair> = {
        let mut out = Vec::new();
        // A smooth junction joins its runs on their shared section; only a
        // corner that turns asks the generators where the walls meet.
        for pair in runs.windows(2) {
            out.push(CornerPair {
                before: pair[0],
                after: pair[1],
                curved: corners.contains(&pair[0].1)
                    && (!straight(pair[0].0, pair[0].1) || !straight(pair[1].0, pair[1].1)),
            });
        }
        if ring && runs.len() > 1 {
            let (before, after) = (runs[runs.len() - 1], runs[0]);
            let turns = stations[before.1]
                .tangent
                .cross(stations[after.0].tangent)
                .magnitude()
                > tol.angular()
                || stations[before.1].tangent.dot(stations[after.0].tangent) < 0.0;
            out.push(CornerPair {
                before,
                after,
                curved: turns && (!straight(before.0, before.1) || !straight(after.0, after.1)),
            });
        }
        out
    };
    let curves: Vec<(ogeom_geom::Curve, (f64, f64), bool)> = {
        let edges: Vec<Shape> = match model.kind_of(spine)? {
            ShapeType::Edge => vec![spine.clone()],
            _ => model.ordered_children_of(spine)?,
        };
        let mut out = Vec::with_capacity(edges.len());
        for edge in &edges {
            let (curve, range) = spine_curve_of(model, edge)?;
            out.push((
                curve,
                range,
                edge.orientation() == ogeom_topo::Orientation::Reversed,
            ));
        }
        out
    };
    let walk = SpineWalk {
        curves,
        stations: &stations,
        normals: &normals,
    };
    let join_reach = tolerance.max(tol.confusion() * 100.0);
    let curved_join = |pair: &CornerPair, ab: (f64, f64)| -> OgeomResult<CornerJoin> {
        let join = walk.join(pair.before, pair.after, ab, tol)?;
        if join.gap > join_reach {
            ogeom_bail!(
                Construction,
                "a skew corner against a curved leg is still owed its frame \
                 law: the legs' generators miss by {}; see docs/PARITY.md, \
                 offset.sweeps",
                join.gap
            );
        }
        Ok(join)
    };
    // Across a curved corner the two legs' skins share the join row as one
    // edge: a run adopts the previous run's end row at its start, and the
    // last run of a cornered ring adopts the first run's start row at its
    // end.
    let shares_start = |ri: usize| -> bool {
        ri > 0
            && corner_pairs
                .iter()
                .any(|pair| pair.curved && pair.after == runs[ri])
    };
    let shares_end = |ri: usize| -> bool {
        ring && ri + 1 == runs.len()
            && corner_pairs
                .iter()
                .any(|pair| pair.curved && pair.before == runs[ri] && pair.after == runs[0])
    };
    // The curved corner a station is a twin of, if any.
    let curved_at = |i: usize| -> Option<&CornerPair> {
        corner_pairs
            .iter()
            .find(|pair| pair.curved && (pair.before.1 == i || pair.after.0 == i))
    };

    // The profile's loops: a face contributes every wire, holes included;
    // a bare wire is one loop.
    let loops: Vec<Shape> = match model.kind_of(profile)? {
        ShapeType::Face => explore(model, profile, Filter::OfType(ShapeType::Wire))?,
        ShapeType::Wire => vec![profile.clone()],
        other => ogeom_bail!(
            Construction,
            "a pipe shell sweeps a planar wire or face, not a {other:?}"
        ),
    };
    if loops.is_empty() {
        ogeom_bail!(Construction, "the profile has no loop to sweep");
    }
    let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
        ogeom_bail!(Construction, "a pipe shell sweeps a planar profile");
    };
    let t0 = stations[0].tangent;
    // Square to the spine's own start tangent, read exactly off its first
    // edge, to an angle's tolerance.
    let exact_t0 = {
        let first = match model.kind_of(spine)? {
            ShapeType::Edge => spine.clone(),
            _ => model.ordered_children_of(spine)?[0].clone(),
        };
        let (curve, range) = spine_curve_of(model, &first)?;
        let reversed = first.orientation() == ogeom_topo::Orientation::Reversed;
        let d = curve.d1_at(if reversed { range.1 } else { range.0 }, tol)?;
        let d = if reversed { -d } else { d };
        d / d.magnitude()
    };
    if plane.normal().vector().cross(exact_t0).magnitude() > tol.angular() {
        ogeom_bail!(
            Construction,
            "the profile leans along its spine; a pipe shell runs square to \
             the start"
        );
    }
    if plane.distance_to(stations[0].at) > tol.confusion() * 100.0 {
        ogeom_bail!(
            Construction,
            "the profile does not sit at the spine's start"
        );
    }

    // Transport: each loop expressed in the start frame's own 2D
    // coordinates, then re-expressed in every station's frame. A loop of one
    // smooth closed edge skins as one wall; a faceted loop skins one strip
    // per edge, cornered at shared vertices, because no single fit can
    // speak a corner.
    let x0 = normals[0];
    let y0 = t0.cross(x0);
    let origin = stations[0].at;
    let flat = |p: Point| -> (f64, f64) { ((p - origin).dot(x0), (p - origin).dot(y0)) };
    // A skew corner against a curved leg: the legs' generators miss, so
    // no join row closes the walls. Such a corner is a mitre instead: the
    // spine is split there, each side swept on straight past the corner
    // and trimmed by the mitre plane, and the pieces fused.
    let skew: Vec<(usize, usize)> = {
        let mut probes: Vec<(f64, f64)> = Vec::new();
        for wire in &loops {
            probes.extend(sample_wire(model, wire, AROUND, tol)?.into_iter().map(flat));
        }
        let mut out = Vec::new();
        for pair in corner_pairs.iter().filter(|pair| pair.curved) {
            let mut worst = 0.0_f64;
            for ab in &probes {
                worst = worst.max(walk.join(pair.before, pair.after, *ab, tol)?.gap);
            }
            if worst > join_reach {
                out.push((pair.before.1, pair.after.0));
            }
        }
        out
    };
    if !skew.is_empty() {
        let probes: Vec<(f64, f64)> = {
            let mut out = Vec::new();
            for wire in &loops {
                out.extend(sample_wire(model, wire, AROUND, tol)?.into_iter().map(flat));
            }
            out
        };
        return mitred_pieces(
            model, profile, spine, &stations, &skew, ring, &probes, tolerance, tol,
        );
    }
    let place = |i: usize, (a, b): (f64, f64)| -> OgeomResult<Point> {
        if let Some(pair) = curved_at(i) {
            return Ok(curved_join(pair, (a, b))?.at);
        }
        let x = normals[i];
        let y = stations[i].tangent.cross(x);
        let p = stations[i].at + x * a + y * b;
        Ok(match mitre[i] {
            Some((corner, n)) => {
                let t = stations[i].tangent;
                p + t * ((corner - p).dot(n) / t.dot(n))
            }
            None => p,
        })
    };
    // The rows a run's skin interpolates. A straight run is its two end
    // rings, ruled: the trimmed prism itself, whether an end is sheared
    // onto a mitre plane or stands on a curved corner's crossing. A curved
    // run is its stations, except that a stretch at a curved corner is
    // re-rowed: each row runs along the generators from the last plain
    // station to the crossing, so the skin is the leg's own surface up to
    // the join and nothing past it.
    let run_rows =
        |(rs, re): (usize, usize), flat_row: &[(f64, f64)]| -> OgeomResult<Vec<Vec<Point>>> {
            let run = (rs, re);
            let joins_at = |at_start: bool| -> OgeomResult<Option<Vec<CornerJoin>>> {
                let pair = corner_pairs.iter().find(|pair| {
                    pair.curved
                        && if at_start {
                            pair.after == run
                        } else {
                            pair.before == run
                        }
                });
                let Some(pair) = pair else {
                    return Ok(None);
                };
                let mut out = Vec::with_capacity(flat_row.len());
                for ab in flat_row {
                    out.push(curved_join(pair, *ab)?);
                }
                Ok(Some(out))
            };
            let (start, end) = (joins_at(true)?, joins_at(false)?);
            #[allow(clippy::cast_precision_loss)]
            let (rsf, ref_) = (rs as f64, re as f64);
            let short = || {
                ogeom_bail!(
                    Construction,
                    "a leg is shorter than its corner's reach; the mitre would \
                 run off its far end"
                )
            };
            let mut rows: Vec<Vec<Point>> = Vec::new();
            if straight(rs, re) {
                if start
                    .as_ref()
                    .is_some_and(|js| js.iter().any(|j| j.s2 >= ref_ - 0.5))
                    || end
                        .as_ref()
                        .is_some_and(|js| js.iter().any(|j| j.s1 <= rsf + 0.5))
                {
                    return short();
                }
                for (at_start, joins) in [(true, &start), (false, &end)] {
                    rows.push(match joins {
                        Some(js) => js.iter().map(|j| j.at).collect(),
                        None => {
                            let i = if at_start { rs } else { re };
                            flat_row
                                .iter()
                                .map(|ab| place(i, *ab))
                                .collect::<OgeomResult<Vec<Point>>>()?
                        }
                    });
                }
                return Ok(rows);
            }
            if start.is_none() && end.is_none() {
                for i in rs..=re {
                    rows.push(
                        flat_row
                            .iter()
                            .map(|ab| place(i, *ab))
                            .collect::<OgeomResult<Vec<Point>>>()?,
                    );
                }
                return Ok(rows);
            }
            // A curved run with a crossing at either end is re-rowed whole:
            // every column runs its own generator from its start to its end,
            // sampled at the same fractions, so the grid's shared parameter is
            // honest for every column; a stretch skewed only near the corner
            // would pace each column differently and the fit would fight it.
            if start
                .as_ref()
                .is_some_and(|js| js.iter().any(|j| j.s2 >= ref_ - 0.5))
                || end
                    .as_ref()
                    .is_some_and(|js| js.iter().any(|j| j.s1 <= rsf + 0.5))
            {
                return short();
            }
            // Each column at equal fractions of its own arc length: the grid's
            // parameter is one for all columns, and a column's pace differs
            // between the leg's curve and its straight extension, so the
            // fractions are taken along the generator, not along its parameter.
            let steps = re - rs;
            let fine = steps * 8;
            let mut columns: Vec<Vec<Point>> = Vec::with_capacity(flat_row.len());
            for (k, ab) in flat_row.iter().enumerate() {
                let lo = start.as_ref().map_or(rsf, |js| js[k].s2);
                let hi = end.as_ref().map_or(ref_, |js| js[k].s1);
                let mut along: Vec<(f64, f64)> = Vec::with_capacity(fine + 1);
                let mut prev: Option<Point> = None;
                let mut length = 0.0;
                for i in 0..=fine {
                    #[allow(clippy::cast_precision_loss)]
                    let sp = lo + (hi - lo) * (i as f64) / (fine as f64);
                    let p = walk.generator(sp, run, *ab, tol)?;
                    if let Some(q) = prev {
                        length += q.distance(p);
                    }
                    along.push((length, sp));
                    prev = Some(p);
                }
                let mut column = Vec::with_capacity(steps + 1);
                for i in 0..=steps {
                    #[allow(clippy::cast_precision_loss)]
                    let target = length * (i as f64) / (steps as f64);
                    let at = along.partition_point(|(l, _)| *l < target).clamp(1, fine);
                    let ((l0, s0), (l1, s1)) = (along[at - 1], along[at]);
                    let f = if l1 > l0 {
                        ((target - l0) / (l1 - l0)).clamp(0.0, 1.0)
                    } else {
                        0.0
                    };
                    column.push(walk.generator(s0 + (s1 - s0) * f, run, *ab, tol)?);
                }
                columns.push(column);
            }
            for i in 0..=steps {
                rows.push(columns.iter().map(|c| c[i]).collect());
            }
            Ok(rows)
        };
    let last = stations.len() - 1;

    enum LoopWall {
        Ring {
            ring0: Shape,
            ring1: Shape,
        },
        Chain {
            bottoms: Vec<Shape>,
            tops: Vec<Shape>,
        },
    }
    let mut faces: Vec<Shape> = Vec::new();
    let mut ends: Vec<LoopWall> = Vec::with_capacity(loops.len());
    for (li, wire) in loops.iter().enumerate() {
        let hole = li != 0;
        let edges = model.ordered_children_of(wire)?;
        let single_smooth = edges.len() == 1 && {
            let (curve, _) = spine_curve_of(model, &edges[0])?;
            !matches!(curve, ogeom_geom::Curve::Line(_))
                && ogeom_algo::edge_vertices(model, &edges[0])?.is_some_and(|(a, b)| a.is_same(&b))
        };
        if single_smooth {
            let samples = sample_wire(model, wire, AROUND, tol)?;
            let flat_row: Vec<(f64, f64)> = samples.iter().map(|p| flat(*p)).collect();
            // One wall per smooth run: a fit across a corner speaks nothing,
            // and the twin stations put both runs' boundary rows on the one
            // mitred ring, where the sew joins them.
            let mut ring0: Option<Shape> = None;
            let mut ring1: Option<Shape> = None;
            for (ri, &(rs, re)) in runs.iter().enumerate() {
                let rows = run_rows((rs, re), &flat_row)?;
                let shared_start = shares_start(ri).then_some(()).and(ring1.as_ref());
                let shared_end = shares_end(ri).then_some(()).and(ring0.as_ref());
                let wall = skinned_wall(model, &rows, (shared_start, shared_end), tolerance, tol)?;
                faces.push(if hole {
                    wall.face.reversed()
                } else {
                    wall.face.clone()
                });
                if ring0.is_none() {
                    ring0 = Some(wall.ring0);
                }
                ring1 = Some(wall.ring1);
            }
            let (Some(ring0), Some(ring1)) = (ring0, ring1) else {
                ogeom_bail!(Construction, "the sweep produced no wall");
            };
            ends.push(LoopWall::Ring { ring0, ring1 });
        } else {
            // Shared corner vertices at both ends of every edge junction.
            let count = edges.len();
            let mut corner_flat: Vec<(f64, f64)> = Vec::with_capacity(count);
            for edge in &edges {
                // Already in the ring's own sense: a reversed edge's
                // vertices come back end first.
                let Some((start, _)) = ogeom_algo::edge_vertices(model, edge)? else {
                    ogeom_bail!(Construction, "a profile edge has no vertices");
                };
                let Some(data) = model.node(&start).and_then(|n| n.data().as_vertex()) else {
                    ogeom_bail!(Construction, "a profile vertex holds no data");
                };
                corner_flat.push(flat(data.point));
            }
            let make_corners = |model: &mut Model, station: usize| -> OgeomResult<Vec<Shape>> {
                let mut out = Vec::with_capacity(corner_flat.len());
                for ab in &corner_flat {
                    out.push(ogeom_algo::make_vertex(model, place(station, *ab)?).shape);
                }
                Ok(out)
            };
            // Corner vertex sets at every run boundary; a kink's twin
            // stations land on the same mitred points, so both runs take
            // the same vertex objects.
            let mut corners_at: Vec<Option<Vec<Shape>>> = vec![None; stations.len()];
            if ring {
                // The wrap is one corner: both runs take the same vertex
                // objects. Every corner's two sheared sections must land on
                // one ring for the loop to close; a planar ring's do
                // exactly, and a skew ring's (whose parallel-carried frame
                // leaves the far tangent's plane) do not, so the residue
                // is measured and the skew ring refused by name rather
                // than sewn hoping.
                let mut worst = 0.0_f64;
                for ab in &corner_flat {
                    worst = worst.max(place(last, *ab)?.distance(place(0, *ab)?));
                    for &k in &kinks {
                        worst = worst.max(place(k, *ab)?.distance(place(k + 1, *ab)?));
                    }
                }
                if worst > tolerance.max(tol.confusion() * 100.0) {
                    ogeom_bail!(
                        Construction,
                        "a skew-cornered ring's sections do not meet on \
                         their mitres; the out-of-plane corner's frame law \
                         is still owed; see docs/PARITY.md, offset.sweeps"
                    );
                }
                let set = make_corners(model, 0)?;
                if worst > tol.confusion() {
                    for v in &set {
                        model.widen(v, ogeom_core::Tolerance::new(worst * 2.0)?)?;
                    }
                }
                corners_at[0] = Some(set.clone());
                corners_at[last] = Some(set);
            } else {
                corners_at[0] = Some(make_corners(model, 0)?);
                corners_at[last] = Some(make_corners(model, last)?);
            }
            for &k in &kinks {
                let set = make_corners(model, k)?;
                corners_at[k] = Some(set.clone());
                corners_at[k + 1] = Some(set);
            }

            // The loop's own centroid line, for orienting each strip.
            let hint_flat = {
                let mut a = 0.0;
                let mut b = 0.0;
                for (fa, fb) in &corner_flat {
                    a += fa;
                    b += fb;
                }
                #[allow(clippy::cast_precision_loss)]
                let n = count as f64;
                (a / n, b / n)
            };

            let mut bottoms = Vec::with_capacity(count);
            let mut tops = Vec::with_capacity(count);
            // Per run: the first strip's start rail, for the last strip to
            // close the loop on, and the previous strip's end rail, for
            // the next to start from: one edge for both, never two fits.
            let mut run_rails: Vec<(Option<Shape>, Option<Shape>)> = vec![(None, None); runs.len()];
            for (ei, edge) in edges.iter().enumerate() {
                let (curve, range) = spine_curve_of(model, edge)?;
                let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
                const ALONG_EDGE: usize = 8;
                let mut flat_row: Vec<(f64, f64)> = Vec::with_capacity(ALONG_EDGE + 1);
                for k in 0..=ALONG_EDGE {
                    #[allow(clippy::cast_precision_loss)]
                    let f = (k as f64) / (ALONG_EDGE as f64);
                    let t = if reversed {
                        range.1 - (range.1 - range.0) * f
                    } else {
                        range.0 + (range.1 - range.0) * f
                    };
                    flat_row.push(flat(curve.point_at(t, tol)?));
                }
                let next = (ei + 1) % count;
                let mut bottom: Option<Shape> = None;
                let mut top: Option<Shape> = None;
                for (ri, &(rs, re)) in runs.iter().enumerate() {
                    let rows = run_rows((rs, re), &flat_row)?;
                    let shared_start = shares_start(ri).then_some(()).and(top.as_ref());
                    let shared_end = shares_end(ri).then_some(()).and(bottom.as_ref());
                    let mid_i = usize::midpoint(rs, re);
                    let hint = {
                        let x = normals[mid_i];
                        let y = stations[mid_i].tangent.cross(x);
                        stations[mid_i].at + x * hint_flat.0 + y * hint_flat.1
                    };
                    let (Some(from), Some(to)) = (&corners_at[rs], &corners_at[re]) else {
                        ogeom_bail!(Construction, "a run boundary has no corners");
                    };
                    let (first_rail0, prev_rail1) = run_rails[ri].clone();
                    let shared_rail0 = if ei > 0 { prev_rail1 } else { None };
                    let shared_rail1 = if ei + 1 == count && count > 1 {
                        first_rail0.clone()
                    } else {
                        None
                    };
                    let strip = skinned_strip(
                        model,
                        &rows,
                        (&from[ei], &from[next], &to[ei], &to[next]),
                        [
                            shared_start,
                            shared_end,
                            shared_rail0.as_ref(),
                            shared_rail1.as_ref(),
                        ],
                        hint,
                        hole,
                        false,
                        tolerance,
                        tol,
                    )?;
                    run_rails[ri] = (
                        if ei == 0 {
                            Some(strip.rail0.clone())
                        } else {
                            first_rail0
                        },
                        Some(strip.rail1.clone()),
                    );
                    faces.push(strip.face.clone());
                    if bottom.is_none() {
                        bottom = Some(strip.bottom);
                    }
                    top = Some(strip.top);
                }
                let (Some(bottom), Some(top)) = (bottom, top) else {
                    ogeom_bail!(Construction, "the sweep produced no strip");
                };
                bottoms.push(bottom);
                tops.push(top);
            }
            ends.push(LoopWall::Chain { bottoms, tops });
        }
    }

    // A cap per end: one plane, one wire per loop, each edge's pcurve the
    // exact projection of its control net into the plane's chart. A ring
    // has no ends: its two boundary rings stand on one mitre plane and the
    // sew joins them.
    for end in 0..if ring { 0 } else { 2 } {
        let (at, outward) = if end == 0 {
            (stations[0].at, -stations[0].tangent)
        } else {
            (stations[last].at, stations[last].tangent)
        };
        let cap_plane = Plane::through(at, Direction::new(outward, tol)?);
        let mut loop_edges: Vec<Vec<Shape>> = Vec::with_capacity(ends.len());
        for wall in &ends {
            loop_edges.push(match wall {
                LoopWall::Ring { ring0, ring1 } => {
                    vec![if end == 0 {
                        ring0.clone()
                    } else {
                        ring1.clone()
                    }]
                }
                LoopWall::Chain { bottoms, tops } => {
                    if end == 0 {
                        bottoms.clone()
                    } else {
                        tops.clone()
                    }
                }
            });
        }
        let mut reach = 1.0_f64;
        for edges in &loop_edges {
            for edge in edges {
                let (curve, range) = spine_curve_of(model, edge)?;
                for t in 0..8 {
                    let p =
                        curve.point_at(range.0 + (range.1 - range.0) * f64::from(t) / 8.0, tol)?;
                    reach = reach.max(p.distance(at) * 2.0);
                }
            }
        }
        let cap_surface: SurfaceGeometry =
            PlaneSurface::over(cap_plane, (-reach, reach), (-reach, reach))?.into();
        let mut wires: Vec<Shape> = Vec::with_capacity(loop_edges.len());
        for edges in &loop_edges {
            wires.push(ogeom_algo::make_wire(model, edges, tol)?.shape);
        }
        let face = ogeom_algo::make_face(model, cap_surface.clone(), &wires, tol)?.shape;
        let cap_id = {
            let Some(node) = model.node(&face) else {
                ogeom_bail!(Dangling, "the cap just built is not in this model");
            };
            let ogeom_topo::NodeData::Face(data) = node.data() else {
                ogeom_bail!(Construction, "the cap holds no face data");
            };
            data.surface
        };
        let frame = cap_plane.frame();
        for edges in &loop_edges {
            for edge in edges {
                let (curve, range) = spine_curve_of(model, edge)?;
                let ogeom_geom::Curve::BSpline(bs) = &curve else {
                    ogeom_bail!(Construction, "a swept ring is not a spline");
                };
                // A planar polynomial spline's chart image is the same-degree
                // spline of the projected control points: affine, so exact.
                let control2: Vec<Point2> = bs
                    .control_points()
                    .iter()
                    .map(|w| {
                        let local = frame.to_local(w.point());
                        Point2::new(local.x, local.y)
                    })
                    .collect();
                let pcurve: ogeom_geom::PlanarCurve =
                    ogeom_geom::BSpline2d::new(bs.knots().clone(), control2, tol)?.into();
                ogeom_algo::attach_pcurve(
                    model,
                    edge,
                    pcurve,
                    cap_id,
                    ogeom_topo::Location::identity(),
                    range,
                )?;
            }
        }
        faces.push(face);
    }

    let sewn = sew(model, &faces, tol)?;
    let mut built = if ring {
        // A holed ring sews into one shell per profile loop: the outer
        // bounds the material, each hole a void tunnel. Largest bound
        // first, the voids' faces already turned at build.
        if sewn.shells.is_empty() {
            ogeom_bail!(Construction, "the pipe shell did not close");
        }
        for shell in &sewn.shells {
            if !ogeom_algo::is_shell_closed(model, shell)? {
                ogeom_bail!(Construction, "the pipe shell did not close");
            }
        }
        let mut ordered = sewn.shells.clone();
        let mut sized: Vec<(f64, Shape)> = Vec::with_capacity(ordered.len());
        for shell in ordered.drain(..) {
            let bound = ogeom_algo::shape_bounds(model, &shell, tol)?;
            sized.push((bound.diagonal(), shell));
        }
        sized.sort_by(|a, b| b.0.partial_cmp(&a.0).unwrap_or(core::cmp::Ordering::Equal));
        let shells: Vec<Shape> = sized.into_iter().map(|(_, s)| s).collect();
        make_solid(model, &shells)?
    } else {
        if sewn.shells.len() != 1 || !ogeom_algo::is_shell_closed(model, &sewn.shells[0])? {
            if std::env::var_os("OGEOM_DEBUG_SWEEP").is_some() {
                eprintln!(
                    "SWEEP: {} shells, {} free edges from {} faces",
                    sewn.shells.len(),
                    sewn.free_edges.len(),
                    faces.len()
                );
                for edge in &sewn.free_edges {
                    let (curve, range) = spine_curve_of(model, edge)?;
                    let a = curve.point_at(range.0, tol)?;
                    let b = curve.point_at(range.1, tol)?;
                    let m = curve.point_at(f64::midpoint(range.0, range.1), tol)?;
                    let t = model.tolerance_of(edge)?.map_or(0.0, |t| t.get());
                    eprintln!(
                        "  free ({:.3},{:.3},{:.3}) -> ({:.3},{:.3},{:.3}) via ({:.3},{:.3},{:.3}) tol {t:.2e}",
                        a.x, a.y, a.z, b.x, b.y, b.z, m.x, m.y, m.z
                    );
                }
            }
            ogeom_bail!(Construction, "the pipe shell did not close");
        }
        make_solid(model, std::slice::from_ref(&sewn.shells[0]))?
    };
    built.history.generate(profile, built.shape.clone());
    built.history.generate(spine, built.shape.clone());
    Ok(built)
}

/// An edge's 3D curve and range, cloned out of the model.
/// A pipe shell whose spine turns a skew corner against a curved leg,
/// built as pieces between such corners and fused.
///
/// With the frame reflected across the mitre plane, a straight leg's walls
/// and a curved leg's cut that plane in sections that differ on the inside
/// of the turn: no single join row closes both. Each piece is swept on
/// straight past its corners (the frame carries unchanged along a straight
/// run), trimmed by each corner's mitre plane, and the pieces fused: where
/// their sections on the plane coincide the caps melt, and where they
/// differ the difference stands as a face of the mitre plane. Exact, and
/// the plain mitre wherever the two sections agree.
#[allow(clippy::too_many_arguments, reason = "one construction, all its data")]
fn mitred_pieces(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    stations: &[SpineStation],
    skew: &[(usize, usize)],
    ring: bool,
    probes: &[(f64, f64)],
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    if model.kind_of(profile)? == ShapeType::Wire {
        // A closed planar wire sweeps the walls of the face it bounds: the
        // solid pieces are mitred and fused as for that face, and its end
        // caps (on the planes square to the spine's ends) are taken off.
        if !ogeom_algo::is_wire_closed(model, profile, tol)? {
            ogeom_bail!(
                Construction,
                "a skew corner against a curved leg is mitred by fusing solid \
                 pieces; an open wire bounds no face to sweep round it"
            );
        }
        let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
            ogeom_bail!(Construction, "a pipe shell sweeps a planar profile");
        };
        let reach = 1e4_f64;
        let surface: SurfaceGeometry =
            PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?.into();
        let face = ogeom_algo::make_face(model, surface, std::slice::from_ref(profile), tol)?.shape;
        let face = realized_profile(model, &face, &Transform::IDENTITY, tol)?;
        let solid = mitred_pieces(
            model, &face, spine, stations, skew, ring, probes, tolerance, tol,
        )?
        .shape;
        let ends: Vec<(Point, Vector)> = if ring {
            Vec::new()
        } else {
            vec![
                (stations[0].at, stations[0].tangent),
                (
                    stations[stations.len() - 1].at,
                    stations[stations.len() - 1].tangent,
                ),
            ]
        };
        let mut walls = Vec::new();
        for f in explore(model, &solid, Filter::OfType(ShapeType::Face))? {
            let Some(ogeom_topo::NodeData::Face(data)) = model.node(&f).map(|n| n.data()) else {
                continue;
            };
            let cap = match model.geometry().surface(data.surface) {
                Some(SurfaceGeometry::Plane(p)) => {
                    let placed = p.plane();
                    ends.iter().any(|(at, n)| {
                        placed.normal().vector().cross(*n).magnitude() <= tol.angular()
                            && placed.distance_to(*at) <= tol.confusion() * 100.0
                    })
                }
                _ => false,
            };
            if !cap {
                walls.push(f);
            }
        }
        let sewn = sew(model, &walls, tol)?;
        let shape = match sewn.shells.as_slice() {
            [shell] => shell.clone(),
            _ => ogeom_algo::make_compound(model, &sewn.shells)?.shape,
        };
        let mut history = History::new();
        history.generate(profile, shape.clone());
        history.generate(spine, shape.clone());
        return Ok(Built::new(shape, history));
    }
    if model.kind_of(profile)? != ShapeType::Face {
        ogeom_bail!(Construction, "a pipe shell sweeps a planar wire or face");
    }
    let edges: Vec<Shape> = match model.kind_of(spine)? {
        ShapeType::Edge => vec![spine.clone()],
        _ => model.ordered_children_of(spine)?,
    };
    let normals = rmf_normals(stations);
    let reach_out = probes
        .iter()
        .map(|(a, b)| a.hypot(*b))
        .fold(0.0_f64, f64::max);

    // Each split corner: the junction's edges, its point, the tangents
    // either side, the frame the far side starts in, and how far each side
    // runs on past it to cover the mitre plane across the whole profile.
    struct Split {
        edge_before: usize,
        at: Point,
        before: Vector,
        after: Vector,
        frame_after: Vector,
        run_on: f64,
    }
    let mut splits: Vec<Split> = Vec::with_capacity(skew.len());
    for &(k, next) in skew {
        let (before, after) = (stations[k].tangent, stations[next].tangent);
        let turn = before.dot(after).clamp(-1.0, 1.0).acos();
        let half = (turn * 0.5).cos();
        if half < 0.05 {
            ogeom_bail!(
                Construction,
                "the spine all but doubles back at a corner; a mitre there \
                 runs off to infinity"
            );
        }
        splits.push(Split {
            edge_before: stations[k].edge,
            at: stations[k].at,
            before,
            after,
            frame_after: normals[next],
            // The mitre plane stands at most `R·tan(φ/2)` past the corner
            // along either leg for a profile reaching `R` from the spine;
            // half as far again clears it with room.
            run_on: reach_out * ((turn * 0.5).tan() * 1.5 + 0.1),
        });
    }
    splits.sort_by_key(|s| s.edge_before);

    // Pieces as runs of spine edges, each between split corners (or the
    // open spine's own ends).
    let count = edges.len();
    let mut pieces: Vec<(Vec<usize>, Option<usize>, Option<usize>)> = Vec::new();
    if ring {
        for (i, split) in splits.iter().enumerate() {
            let next = &splits[(i + 1) % splits.len()];
            let mut run = Vec::new();
            let mut e = (split.edge_before + 1) % count;
            loop {
                run.push(e);
                if e == next.edge_before {
                    break;
                }
                e = (e + 1) % count;
            }
            pieces.push((run, Some(i), Some((i + 1) % splits.len())));
        }
    } else {
        let mut first = 0;
        for (i, split) in splits.iter().enumerate() {
            pieces.push((
                (first..=split.edge_before).collect(),
                i.checked_sub(1),
                Some(i),
            ));
            first = split.edge_before + 1;
        }
        pieces.push(((first..count).collect(), splits.len().checked_sub(1), None));
    }

    let x0 = normals[0];
    let start_frame = Frame::new(
        stations[0].at,
        Direction::new(stations[0].tangent, tol)?,
        Direction::new(x0, tol)?,
        tol,
    )?;
    // One block per split, standing on the mitre plane on the far side of
    // the corner: the piece before the corner is cut by it and the piece
    // after keeps what it shares with it, so both sides' caps are pieces
    // of the block's one face, on one surface and one chart, which is what
    // lets the fuse melt them.
    let mut blocks: Vec<Shape> = Vec::with_capacity(splits.len());
    for split in &splits {
        let n = (split.before + split.after) / (split.before + split.after).magnitude();
        let normal = Direction::new(n, tol)?;
        let plane = Plane::through(split.at, normal);
        let reach = (split.run_on + reach_out) * 4.0;
        let frame = plane.frame();
        let (u, v) = (frame.x().vector(), frame.y().vector());
        let corners: Vec<Point> = [(-1.0, -1.0), (1.0, -1.0), (1.0, 1.0), (-1.0, 1.0)]
            .iter()
            .map(|(a, b)| split.at + u * (a * reach) + v * (b * reach))
            .collect();
        let wire = ogeom_algo::make_polygon(model, &corners, true, tol)?.shape;
        let edges = explore(model, &wire, Filter::OfType(ShapeType::Edge))?;
        let surface: SurfaceGeometry = PlaneSurface::over(
            plane,
            (-reach * 2.0, reach * 2.0),
            (-reach * 2.0, reach * 2.0),
        )?
        .into();
        let base = ogeom_algo::make_face_with_pcurves(model, surface, &[edges], tol)?.shape;
        let block = ogeom_algo::make_prism(model, &base, n * (reach * 2.0), tol)?.shape;
        blocks.push(block);
    }

    let mut joined: Vec<Shape> = Vec::new();
    for (run, start, end) in pieces {
        let mut wire_edges: Vec<Shape> = Vec::new();
        let traversal = |model: &Model, e: usize, at_start: bool| -> OgeomResult<Shape> {
            // In the spine's own sense: a reversed edge's vertices come
            // back end first.
            let Some((a, b)) = ogeom_algo::edge_vertices(model, &edges[e])? else {
                ogeom_bail!(Construction, "a spine edge has no vertices");
            };
            Ok(if at_start { a } else { b })
        };
        if let Some(i) = start {
            let split = &splits[i];
            let far = ogeom_algo::make_vertex(model, split.at - split.after * split.run_on).shape;
            let near = traversal(model, run[0], true)?;
            let line: ogeom_geom::Curve =
                LineCurve::segment(split.at - split.after * split.run_on, split.at, tol)?.into();
            let domain = line.domain();
            wire_edges
                .push(ogeom_algo::make_edge_between(model, line, domain, &far, &near, tol)?.shape);
        }
        wire_edges.extend(run.iter().map(|&e| edges[e].clone()));
        if let Some(i) = end {
            let split = &splits[i];
            let near = traversal(model, run[run.len() - 1], false)?;
            let far = ogeom_algo::make_vertex(model, split.at + split.before * split.run_on).shape;
            let line: ogeom_geom::Curve =
                LineCurve::segment(split.at, split.at + split.before * split.run_on, tol)?.into();
            let domain = line.domain();
            wire_edges
                .push(ogeom_algo::make_edge_between(model, line, domain, &near, &far, tol)?.shape);
        }
        let sub_spine = ogeom_algo::make_wire(model, &wire_edges, tol)?.shape;
        // The profile where this piece starts: the spine's own start keeps
        // the caller's; a piece starting past a corner takes the profile
        // moved into the frame the corner's far side starts in, set back
        // along its run-on.
        let placed = match start {
            None => profile.clone(),
            Some(i) => {
                let split = &splits[i];
                let target = Frame::new(
                    split.at - split.after * split.run_on,
                    Direction::new(split.after, tol)?,
                    Direction::new(split.frame_after, tol)?,
                    tol,
                )?;
                let motion = Transform::from_frame(&target) * Transform::to_frame(&start_frame);
                realized_profile(model, profile, &motion, tol)?
            }
        };
        let mut piece = make_pipe_shell(model, &placed, &sub_spine, false, tolerance, tol)?.shape;
        if model.kind_of(&piece)? != ShapeType::Solid {
            ogeom_bail!(Construction, "a mitred piece did not sweep into a solid");
        }
        if let Some(i) = start {
            piece = ogeom_bool::common(model, &piece, &blocks[i], tol)?.shape;
        }
        if let Some(i) = end {
            piece = ogeom_bool::cut(model, &piece, &blocks[i], tol)?.shape;
        }
        joined.push(piece);
    }
    let result = fuse_in_order(model, joined, tol)?;
    let Some(shape) = result else {
        ogeom_bail!(Construction, "the spine produced no piece to sweep");
    };
    let mut history = History::new();
    history.generate(profile, shape.clone());
    for edge in &edges {
        history.generate(edge, shape.clone());
    }
    Ok(Built::new(shape, history))
}

/// A planar profile face rebuilt under a rigid motion: every edge's curve
/// moved and re-bounded, vertices shared, the face on the moved plane.
fn realized_profile(
    model: &mut Model,
    profile: &Shape,
    motion: &Transform,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    realized_profile_wound(model, profile, motion, None, tol)
}

/// As [`realized_profile`], each ring wound about `about` where given: the
/// outer ring turning positively, every hole the other way.
fn realized_profile_wound(
    model: &mut Model,
    profile: &Shape,
    motion: &Transform,
    about: Option<Vector>,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    use ogeom_geom::Transformable as _;
    let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
        ogeom_bail!(Construction, "a pipe shell sweeps a planar profile");
    };
    let moved_plane = Plane::through(
        motion.apply(plane.origin()),
        Direction::new(motion.apply_vector(plane.normal().vector()), tol)?,
    );
    let mut vertices: std::collections::HashMap<ogeom_topo::TShapeId, Shape> =
        std::collections::HashMap::new();
    let mut edge_copies: std::collections::HashMap<ogeom_topo::TShapeId, Shape> =
        std::collections::HashMap::new();
    let mut wires: Vec<Vec<Shape>> = Vec::new();
    for wire in explore(model, profile, Filter::OfType(ShapeType::Wire))? {
        let mut ring = Vec::new();
        for edge in model.ordered_children_of(&wire)? {
            let copy = match edge_copies.get(&edge.node()) {
                Some(done) => done.clone(),
                None => {
                    let (curve, range) = spine_curve_of(model, &edge)?;
                    let placed = curve.transformed(&edge.transform(model.datums())?, tol)?;
                    let moved = placed.transformed(motion, tol)?;
                    // The copy is of the edge itself, in its own sense; the
                    // ring's use of it is reapplied below.
                    let own = if edge.orientation() == ogeom_topo::Orientation::Reversed {
                        edge.reversed()
                    } else {
                        edge.clone()
                    };
                    let Some((a, b)) = ogeom_algo::edge_vertices(model, &own)? else {
                        ogeom_bail!(Construction, "a profile edge has no vertices");
                    };
                    let mut ends = Vec::with_capacity(2);
                    for v in [a, b] {
                        let key = v.node();
                        let held = match vertices.get(&key) {
                            Some(done) => done.clone(),
                            None => {
                                let Some(data) = model.node(&v).and_then(|n| n.data().as_vertex())
                                else {
                                    ogeom_bail!(Construction, "a profile vertex holds no data");
                                };
                                let at = v.transform(model.datums())?.apply(data.point);
                                let fresh = ogeom_algo::make_vertex(model, motion.apply(at)).shape;
                                vertices.insert(key, fresh.clone());
                                fresh
                            }
                        };
                        ends.push(held);
                    }
                    let fresh = ogeom_algo::make_edge_between(
                        model, moved, range, &ends[0], &ends[1], tol,
                    )?
                    .shape;
                    edge_copies.insert(edge.node(), fresh.clone());
                    fresh
                }
            };
            ring.push(if edge.orientation() == ogeom_topo::Orientation::Reversed {
                copy.reversed()
            } else {
                copy
            });
        }
        if let Some(axis) = about {
            let turning = ring_turning(model, &ring, axis, tol)?;
            let outer = wires.is_empty();
            if (turning > 0.0) != outer {
                ring = ring.iter().rev().map(Shape::reversed).collect();
            }
        }
        wires.push(ring);
    }
    let reach = 1e4_f64;
    let surface: SurfaceGeometry =
        PlaneSurface::over(moved_plane, (-reach, reach), (-reach, reach))?.into();
    Ok(ogeom_algo::make_face_with_pcurves(model, surface, &wires, tol)?.shape)
}

/// Twice the signed area a ring of edges encloses about `axis`, from its
/// edges sampled in the ring's own sense.
fn ring_turning(model: &Model, ring: &[Shape], axis: Vector, tol: Tolerances) -> OgeomResult<f64> {
    let mut points: Vec<Point> = Vec::new();
    for edge in ring {
        let (curve, range) = spine_curve_of(model, edge)?;
        let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
        for i in 0..32 {
            let f = f64::from(i) / 32.0;
            let t = if reversed {
                range.1 - (range.1 - range.0) * f
            } else {
                range.0 + (range.1 - range.0) * f
            };
            points.push(curve.point_at(t, tol)?);
        }
    }
    let Some(&origin) = points.first() else {
        return Ok(0.0);
    };
    let n = points.len();
    Ok((0..n)
        .map(|i| {
            (points[i] - origin)
                .cross(points[(i + 1) % n] - origin)
                .dot(axis)
        })
        .sum())
}

pub(crate) fn spine_curve_of(
    model: &Model,
    edge: &Shape,
) -> OgeomResult<(ogeom_geom::Curve, (f64, f64))> {
    let Some(data) = model.node(edge).and_then(|n| n.data().as_edge()) else {
        ogeom_bail!(Construction, "an edge holds no data");
    };
    let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
        ogeom_bail!(Construction, "an edge has no curve");
    };
    let Some(geometry) = model.geometry().curve(*curve) else {
        ogeom_bail!(Dangling, "curve is not in this model");
    };
    Ok((geometry.clone(), *range))
}

/// The pipe shell around a spine that loops back on itself: one wall,
/// closed both ways round, no caps at all.
///
/// The frames are rotation-minimizing with the loop's holonomy paid off:
/// transported round a closed spine, the frame comes home twisted by some
/// angle, and that twist is spread back along the arc so the last station's
/// frame *is* the first's; without it the closed fit fights a helical
/// grid. The profile may be smooth or faceted, and each hole sweeps a void
/// tunnel of its own shell; the Frenet law rides the loop as well.
fn closed_pipe_shell(
    model: &mut Model,
    profile: &Shape,
    spine: &Shape,
    mut stations: Vec<SpineStation>,
    frenet: bool,
    tolerance: f64,
    tol: Tolerances,
) -> OgeomResult<Built> {
    // The walk visits the join twice; the loop owns it once.
    stations.pop();
    if stations.len() < 3 {
        ogeom_bail!(Construction, "a closed spine needs room to turn");
    }
    // A sharp corner is a kink, and a kinked ring is mitred by the caller;
    // one arriving here has a heading that jumps between two stations.
    for i in 0..stations.len() {
        let next = &stations[(i + 1) % stations.len()];
        if stations[i].tangent.dot(next.tangent) < 0.9 {
            ogeom_bail!(
                Construction,
                "a closed spine turns too sharply between two of its stations to skin"
            );
        }
    }

    let loops: Vec<Shape> = match model.kind_of(profile)? {
        ShapeType::Face => explore(model, profile, Filter::OfType(ShapeType::Wire))?,
        ShapeType::Wire => vec![profile.clone()],
        other => ogeom_bail!(
            Construction,
            "a pipe shell sweeps a planar wire or face, not a {other:?}"
        ),
    };
    // Every wire sweeps its own closed shell: the outer boundary first,
    // each hole a void tunnel inside it.
    let profile_loop = &loops[0];
    let edges = model.ordered_children_of(profile_loop)?;
    let smooth = edges.len() == 1
        && ogeom_algo::edge_vertices(model, &edges[0])?.is_some_and(|(a, b)| a.is_same(&b));
    let Some(plane) = ogeom_algo::find_plane(model, profile, tol)? else {
        ogeom_bail!(Construction, "a pipe shell sweeps a planar profile");
    };
    let t0 = stations[0].tangent;
    // Square to the spine's own start tangent, read exactly off its first
    // edge, to an angle's tolerance.
    let exact_t0 = {
        let first = match model.kind_of(spine)? {
            ShapeType::Edge => spine.clone(),
            _ => model.ordered_children_of(spine)?[0].clone(),
        };
        let (curve, range) = spine_curve_of(model, &first)?;
        let reversed = first.orientation() == ogeom_topo::Orientation::Reversed;
        let d = curve.d1_at(if reversed { range.1 } else { range.0 }, tol)?;
        let d = if reversed { -d } else { d };
        d / d.magnitude()
    };
    if plane.normal().vector().cross(exact_t0).magnitude() > tol.angular() {
        ogeom_bail!(
            Construction,
            "the profile leans along its spine; a pipe shell runs square to \
             the start"
        );
    }
    if plane.distance_to(stations[0].at) > tol.confusion() * 100.0 {
        ogeom_bail!(
            Construction,
            "the profile does not sit at the spine's start"
        );
    }

    // Frames with the loop's mismatch paid off: carry once more back to the
    // start, read the twist between departure and return, and spread it
    // along the arc. Rotation-minimizing frames owe this for their
    // holonomy; the Frenet law owes it too, because straight stretches
    // carry the frame through by continuation and the continuation is
    // path-dependent. One reconciliation serves both.
    let mut normals: Vec<Vector> = if frenet {
        // The Frenet frame is the spine's own, single-valued round a loop:
        // read with wrapped neighbours it closes on itself and owes no
        // reconciliation. Read from a walk that visits the join twice it
        // does not: the one-sided differences at the walk's two ends
        // disagree with the interior, and the strips built on them miss
        // each other at the join by that kink.
        frenet_normals_closed(&stations, tol)?
    } else {
        let mut extended = stations.clone();
        extended.push(stations[0]);
        let carried = rmf_normals(&extended);
        let (n0, n_home) = (carried[0], carried[carried.len() - 1]);
        let twist = (n0.cross(n_home).dot(t0)).atan2(n0.dot(n_home));
        let mut lengths = vec![0.0_f64];
        for pair in extended.windows(2) {
            let last = lengths[lengths.len() - 1];
            lengths.push(last + pair[0].at.distance(pair[1].at));
        }
        let total = lengths[lengths.len() - 1];
        carried
            .iter()
            .take(stations.len())
            .enumerate()
            .map(|(i, n)| {
                let phi = -twist * lengths[i] / total;
                let t = extended[i].tangent;
                *n * phi.cos() + t.cross(*n) * phi.sin()
            })
            .collect()
    };
    for (n, station) in normals.iter_mut().zip(&stations) {
        // Re-square each corrected normal against its own tangent.
        let v = *n - station.tangent * n.dot(station.tangent);
        *n = v / v.magnitude();
    }

    let x0 = normals[0];
    let origin = stations[0].at;
    let mut shells: Vec<Shape> = Vec::with_capacity(loops.len());
    for (li, wire) in loops.iter().enumerate() {
        let wire_edges = model.ordered_children_of(wire)?;
        let wire_smooth = wire_edges.len() == 1
            && ogeom_algo::edge_vertices(model, &wire_edges[0])?
                .is_some_and(|(a, b)| a.is_same(&b));
        let shell = closed_loop_shell(
            model,
            wire,
            &wire_edges,
            wire_smooth,
            &stations,
            &normals,
            (origin, x0),
            tolerance,
            tol,
        )?;
        // The material side follows each loop's own winding, so it is
        // read off the shell itself: the outer shell faces out of what it
        // encloses, a void's faces toward its tunnel.
        let enclosed = shell_signed_volume(model, &shell, tol)?;
        let outward = enclosed > 0.0;
        shells.push(if outward == (li == 0) {
            shell
        } else {
            shell.reversed()
        });
    }
    let mut built = make_solid(model, &shells)?;
    built.history.generate(profile, built.shape.clone());
    built.history.generate(spine, built.shape.clone());
    let _ = (smooth, profile_loop, edges);
    Ok(built)
}

/// The volume a closed shell encloses as it faces, from its mesh: negative
/// where its faces point into what it bounds.
fn shell_signed_volume(model: &Model, shell: &Shape, tol: Tolerances) -> OgeomResult<f64> {
    let mesh = ogeom_mesh::triangulate(model, shell, ogeom_mesh::Deflection::default(), tol)?;
    Ok(mesh
        .triangles
        .iter()
        .map(|t| {
            let [a, b, c] = t.map(|i| mesh.positions[i as usize].to_vector());
            a.dot(b.cross(c)) / 6.0
        })
        .sum())
}

/// Sample a spine (one edge or a wire of them) into stations, each edge
/// given a station count by its own turning.
fn shell_stations(model: &Model, spine: &Shape, tol: Tolerances) -> OgeomResult<Vec<SpineStation>> {
    let edges: Vec<Shape> = match model.kind_of(spine)? {
        ShapeType::Edge => vec![spine.clone()],
        ShapeType::Wire => model.ordered_children_of(spine)?,
        other => ogeom_bail!(
            Construction,
            "a pipe shell runs along an edge or a wire, not a {other:?}"
        ),
    };
    if edges.is_empty() {
        ogeom_bail!(Construction, "the spine has no edge to run along");
    }
    let mut stations: Vec<SpineStation> = Vec::new();
    for (ei, edge) in edges.iter().enumerate() {
        let (curve, range) = {
            let Some(data) = model.node(edge).and_then(|n| n.data().as_edge()) else {
                ogeom_bail!(Construction, "a spine edge holds no data");
            };
            let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
                ogeom_bail!(Construction, "a spine edge has no curve");
            };
            let Some(geometry) = model.geometry().curve(*curve) else {
                ogeom_bail!(Dangling, "curve is not in this model");
            };
            (geometry.clone(), *range)
        };
        let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
        // Stations by turning: sample tangents coarsely, sum the angles, and
        // give each edge enough stations that no step turns more than a few
        // degrees. A straight edge keeps a healthy minimum for the fit.
        let turning = {
            let mut sum = 0.0_f64;
            let mut last: Option<Vector> = None;
            for i in 0..=16 {
                let t = range.0 + (range.1 - range.0) * f64::from(i) / 16.0;
                let d = curve.d1_at(t, tol)?;
                let m = d.magnitude();
                if m <= tol.confusion() {
                    continue;
                }
                let u = d / m;
                if let Some(prev) = last {
                    sum += prev.dot(u).clamp(-1.0, 1.0).acos() * 16.0 / 16.0;
                }
                last = Some(u);
            }
            sum
        };
        #[allow(clippy::cast_possible_truncation, clippy::cast_sign_loss)]
        let count = (turning / (core::f64::consts::TAU / 64.0)).ceil().max(8.0) as usize;
        for i in 0..=count {
            #[allow(clippy::cast_precision_loss)]
            let f = (i as f64) / (count as f64);
            let t = if reversed {
                range.1 - (range.1 - range.0) * f
            } else {
                range.0 + (range.1 - range.0) * f
            };
            let p = curve.point_at(t, tol)?;
            let d = curve.d1_at(t, tol)?;
            let m = d.magnitude();
            if m <= tol.confusion() {
                ogeom_bail!(Construction, "the spine is degenerate at {t}");
            }
            let tangent = if reversed { -(d / m) } else { d / m };
            if let Some(prev) = stations.last()
                && prev.edge == ei
                && prev.at.distance(p) <= tol.confusion()
                && prev.tangent.cross(tangent).magnitude() <= tol.angular()
                && prev.tangent.dot(tangent) > 0.0
            {
                continue;
            }
            // A station coincident with the last but on the next edge is a
            // twin: heading elsewhere it is a *corner* the sweep mitres, and
            // heading on it is a smooth junction where the next edge's own
            // run of skin begins, since one fit across two curves' joins
            // cannot follow the step in their curvature.
            stations.push(SpineStation {
                at: p,
                tangent,
                edge: ei,
                t,
            });
        }
    }
    if stations.len() < 2 {
        ogeom_bail!(Construction, "the spine collapses to a point");
    }
    Ok(stations)
}

/// Frenet normals: each station's frame turns with the spine's own
/// curvature, read from the tangents' finite differences. Straight runs
/// carry the last bending station's normal forward; a spine that never
/// bends has no Frenet frame at all and is refused by name.
fn frenet_normals(stations: &[SpineStation], tol: Tolerances) -> OgeomResult<Vec<Vector>> {
    let mut normals: Vec<Option<Vector>> = Vec::with_capacity(stations.len());
    for i in 0..stations.len() {
        let (before, after) = (
            &stations[i.saturating_sub(1)],
            &stations[(i + 1).min(stations.len() - 1)],
        );
        let dt = after.tangent - before.tangent;
        let t = stations[i].tangent;
        let bend = dt - t * dt.dot(t);
        let m = bend.magnitude();
        normals.push(if m > tol.angular().max(1e-9) {
            Some(bend / m)
        } else {
            None
        });
    }
    // Carry forward, then backward, so straight lead-ins take the first
    // bend's frame rather than none.
    let mut carried: Vec<Vector> = Vec::with_capacity(stations.len());
    let mut last: Option<Vector> = None;
    for n in &normals {
        if let Some(n) = n {
            last = Some(*n);
        }
        carried.push(last.unwrap_or(Vector::new(0.0, 0.0, 0.0)));
    }
    let mut ahead: Option<Vector> = None;
    for i in (0..stations.len()).rev() {
        if let Some(n) = normals[i] {
            ahead = Some(n);
        } else if carried[i].magnitude() < 0.5
            && let Some(n) = ahead
        {
            carried[i] = n;
        }
    }
    if carried.iter().any(|n| n.magnitude() < 0.5) {
        ogeom_bail!(
            Construction,
            "a straight spine has no Frenet frame; use the \
             rotation-minimizing default"
        );
    }
    Ok(carried)
}

/// Frenet normals round a closed loop: each station's bend read from its
/// neighbours across the join as well, so the field is periodic. A loop
/// with a straight stretch carries the last bend's normal through it, as
/// the open form does; a loop that never bends has no Frenet frame.
fn frenet_normals_closed(stations: &[SpineStation], tol: Tolerances) -> OgeomResult<Vec<Vector>> {
    let n = stations.len();
    let mut normals: Vec<Option<Vector>> = Vec::with_capacity(n);
    for i in 0..n {
        let (before, after) = (&stations[(i + n - 1) % n], &stations[(i + 1) % n]);
        let dt = after.tangent - before.tangent;
        let t = stations[i].tangent;
        let bend = dt - t * dt.dot(t);
        let m = bend.magnitude();
        normals.push(if m > tol.angular().max(1e-9) {
            Some(bend / m)
        } else {
            None
        });
    }
    let Some(first_bend) = normals.iter().position(Option::is_some) else {
        ogeom_bail!(
            Construction,
            "a straight spine has no Frenet frame; use the \
             rotation-minimizing default"
        );
    };
    // Carry forward round the loop from the first bend, so a straight
    // stretch anywhere takes the bend behind it.
    let mut carried: Vec<Vector> = vec![Vector::new(0.0, 0.0, 0.0); n];
    let mut last = normals[first_bend].unwrap_or_else(|| unreachable!());
    for k in 0..n {
        let i = (first_bend + k) % n;
        if let Some(bend) = normals[i] {
            last = bend;
        }
        carried[i] = last;
    }
    Ok(carried)
}

// --- the evolved shape -------------------------------------------------------

/// One station of the spine: where it is, which way it runs, and how it gets
/// there.
struct Station {
    /// Where the traversal enters and leaves this edge.
    from: Point,
    to: Point,
    /// The unit tangent at each end, in the direction of travel.
    tangent_in: Vector,
    tangent_out: Vector,
    /// A straight run, or a turn about an axis through an angle.
    turn: Option<(ogeom_math::Axis, f64)>,
}

/// Sweep a profile along a spine, the way a moulding runs round a frame.
///
/// The spine is a **planar** wire, or a planar face whose outer wire is taken.
/// The profile is a wire standing in a plane that contains the spine's own
/// normal, positioned where the spine starts. What comes back is what the
/// profile sweeps out as it travels the spine, always square to it:
///
/// - a straight spine edge extrudes the profile (a prism);
/// - a circular one turns it about that arc's own axis (a revolution);
/// - and each corner between them turns it about the corner, through exactly
///   the angle the spine turns there, which is the join the 2D offset makes
///   for the same reason.
///
/// Every piece is exact: the surfaces are the ones a prism and a revolution
/// give for the profile's own curves, and nothing is fitted. The pieces are
/// then unioned, which is the assembly's real name: consecutive pieces meet
/// on the *same* placed profile, and a coincident face is what the boolean
/// identifies rather than probes across.
///
/// # Volume or shell
///
/// The result is always a volume, and which spine is given is what says
/// whether there is one to have. A **closed** profile bounds its own section
/// and sweeps a solid along either kind of spine. An **open** one does not,
/// and there is exactly one honest way to close it: against the plane a
/// **face** spine was drawn in, whose own plane the profile's two ends must
/// reach. An open profile along a wire spine is refused, and the refusal says
/// which spine would close it.
///
/// # Errors
///
/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the spine is
/// not a planar wire or face, if it carries an edge that is neither straight
/// nor circular, if the profile is not planar, if the profile's plane does not
/// contain the spine's normal or does not cut across it (a profile that leans
/// or lies along is not square to the spine), or if an open profile has no
/// spine plane to close against.
/// [`OgeomError::NotDone`](ogeom_core::OgeomError::NotDone) where a corner's turn
/// would sweep the profile across the corner itself, which no revolution can
/// express.
pub fn make_evolved(
    model: &mut Model,
    spine: &Shape,
    profile: &Shape,
    tol: Tolerances,
) -> OgeomResult<Built> {
    use ogeom_algo::{make_prism, make_revolution, transformed};

    let (wire, capped_by_spine_plane) = match model.kind_of(spine)? {
        ShapeType::Face => {
            let wires = explore(model, spine, Filter::OfType(ShapeType::Wire))?;
            let Some(outer) = wires.first().cloned() else {
                ogeom_bail!(Construction, "a face with no wire has no spine to run");
            };
            (outer, true)
        }
        ShapeType::Wire => (spine.clone(), false),
        other => ogeom_bail!(
            Construction,
            "a {other:?} is not a spine; sweep along a wire or a planar face"
        ),
    };

    let stations = spine_stations(model, &wire, tol)?;
    if stations.is_empty() {
        ogeom_bail!(Construction, "a spine with no edges goes nowhere");
    }
    let normal = spine_normal(&stations, tol)?;
    let (profile_origin, profile_normal) = profile_plane(model, profile, tol)?;
    if profile_normal.dot(normal.vector()).abs() > tol.angular() {
        ogeom_bail!(
            Construction,
            "the profile's plane must contain the spine's normal, or the \
             profile is not square to the spine it travels"
        );
    }
    let start = &stations[0];
    if profile_normal.cross(start.tangent_in).magnitude() > tol.angular() {
        ogeom_bail!(
            Construction,
            "the profile's plane must cut the spine across, not run along it: \
             the profile is not square to the spine it travels"
        );
    }
    let _ = profile_origin;
    let reference = (start.from, start.tangent_in);

    // The profile as a face, which is what makes each swept piece a *solid*
    // and the assembly a union rather than a hopeful sew. An open profile is
    // closed against the spine's own plane, which is exactly what a face
    // spine offers and a wire spine does not.
    let section = profile_face(
        model,
        profile,
        profile_normal,
        capped_by_spine_plane.then(|| Plane::through(start.from, normal)),
        tol,
    )?;

    let mut pieces: Vec<Shape> = Vec::new();
    for (index, station) in stations.iter().enumerate() {
        // The corner *before* this station, so the pieces come out in the
        // order the spine runs them.
        if index > 0 {
            let previous = &stations[index - 1];
            if let Some(piece) = corner_piece(
                model,
                &section,
                reference,
                previous.to,
                previous.tangent_out,
                station.tangent_in,
                normal,
                tol,
            )? {
                pieces.push(piece);
            }
        }
        let placed = transformed(
            model,
            &section,
            station_transform(reference, station.from, station.tangent_in, normal, tol)?,
        )?
        .shape;
        pieces.push(match station.turn {
            None => make_prism(model, &placed, station.to - station.from, tol)?.shape,
            Some((axis, angle)) => make_revolution(model, &placed, axis, angle, tol)?.shape,
        });
    }
    // A closed spine turns at the join between its last edge and its first
    // just as it does anywhere else.
    let last = &stations[stations.len() - 1];
    if last.to.distance(start.from) <= tol.confusion()
        && let Some(piece) = corner_piece(
            model,
            &section,
            reference,
            last.to,
            last.tangent_out,
            start.tangent_in,
            normal,
            tol,
        )?
    {
        pieces.push(piece);
    }

    // The union, in the order the spine runs: consecutive pieces meet on the
    // *same* placed profile, which is the coincident-face case the boolean
    // resolves by identifying it rather than by probing across it.
    let mut history = History::new();
    let Some(shape) = fuse_in_order(model, pieces, tol)? else {
        ogeom_bail!(Construction, "the spine produced no piece to sweep");
    };
    history.generate(spine, shape.clone());
    history.generate(profile, shape.clone());
    Ok(Built::new(shape, history))
}

/// The profile as a face.
///
/// A closed profile bounds its own area. An open one does not, and there is
/// exactly one honest way to close it: against the plane the spine was given
/// as a face *in*, which is what a face spine says to do and a wire spine has
/// no answer for. The closing segment runs between the profile's two ends, and
/// both have to be on that plane or the profile does not reach it.
fn profile_face(
    model: &mut Model,
    profile: &Shape,
    profile_normal: Vector,
    against: Option<Plane>,
    tol: Tolerances,
) -> OgeomResult<Shape> {
    if model.kind_of(profile)? == ShapeType::Face {
        return Ok(profile.clone());
    }
    if model.kind_of(profile)? != ShapeType::Wire {
        ogeom_bail!(Construction, "a profile is a wire or a face");
    }
    let mut edges = ogeom_topo::explore(model, profile, Filter::OfType(ShapeType::Edge))?;
    let closed = ogeom_algo::is_wire_closed(model, profile, tol)?;
    if !closed {
        let Some(plane) = against else {
            ogeom_bail!(
                Construction,
                "an open profile sweeps a shell, not a volume; give the spine \
                 as a planar face for its plane to close the profile against, \
                 or close the profile itself"
            );
        };
        let [(from, v0), (to, v1)] = wire_ends(model, profile, tol)?;
        for end in [from, to] {
            if plane.signed_distance_to(end).abs() > tol.confusion() * 1e2 {
                ogeom_bail!(
                    Construction,
                    "an open profile is closed against the spine face's own \
                     plane, and this one does not reach it"
                );
            }
        }
        // Built on the profile's *own* end vertices, so the closed ring is a
        // wire rather than edges that merely touch.
        let line = LineCurve::new(ogeom_math::Axis {
            location: from,
            direction: Direction::new(to - from, tol)?,
        });
        edges.push(
            make_edge_between(model, line.into(), (0.0, from.distance(to)), &v0, &v1, tol)?.shape,
        );
    }
    let ordered = ogeom_algo::order_edges(model, &edges, tol)?;
    let mut bound = ogeom_math::Aabb::EMPTY;
    for edge in &ordered {
        bound = bound.union(&ogeom_algo::shape_bounds(model, edge, tol)?);
    }
    let Some(centre) = bound.centre() else {
        ogeom_bail!(Construction, "a profile with no extent sweeps nothing");
    };
    let reach = bound.diagonal().mul_add(2.0, 1.0);
    let plane = Plane::through(centre, Direction::new(profile_normal, tol)?);
    let surface = PlaneSurface::over(plane, (-reach, reach), (-reach, reach))?;
    Ok(make_face_with_pcurves(model, surface.into(), &[ordered], tol)?.shape)
}

/// Where an open wire begins and ends: the point, and the vertex there.
fn wire_ends(model: &Model, wire: &Shape, tol: Tolerances) -> OgeomResult<[(Point, Shape); 2]> {
    let mut counts: Vec<(Point, Shape, usize)> = Vec::new();
    for edge in explore(model, wire, Filter::OfType(ShapeType::Edge))? {
        for v in explore(model, &edge, Filter::OfType(ShapeType::Vertex))? {
            let Some(data) = model.node(&v).and_then(|n| n.data().as_vertex()) else {
                continue;
            };
            let at = v.transform(model.datums())?.apply(data.point);
            match counts
                .iter_mut()
                .find(|(p, _, _)| p.distance(at) <= tol.confusion() * 10.0)
            {
                Some((_, _, n)) => *n += 1,
                None => counts.push((at, v.clone(), 1)),
            }
        }
    }
    let free: Vec<(Point, Shape)> = counts
        .into_iter()
        .filter(|(_, _, n)| *n == 1)
        .map(|(p, v, _)| (p, v))
        .collect();
    if free.len() != 2 {
        ogeom_bail!(
            Construction,
            "an open profile has exactly two ends; this one has {}",
            free.len()
        );
    }
    let mut ends = free.into_iter();
    let (Some(a), Some(b)) = (ends.next(), ends.next()) else {
        ogeom_bail!(Construction, "the profile lost an end between checks");
    };
    Ok([a, b])
}

/// The spine, edge by edge, in the order the wire runs it.
fn spine_stations(model: &Model, wire: &Shape, tol: Tolerances) -> OgeomResult<Vec<Station>> {
    let mut out = Vec::new();
    for edge in explore(model, wire, Filter::OfType(ShapeType::Edge))? {
        let Some(data) = model.node(&edge).and_then(|n| n.data().as_edge()) else {
            ogeom_bail!(Construction, "a spine edge is not in this model");
        };
        let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
            ogeom_bail!(Construction, "a spine edge with no curve runs nowhere");
        };
        let Some(geometry) = model.geometry().curve(*curve) else {
            ogeom_bail!(Dangling, "curve is not in this model");
        };
        let placed = geometry.transformed(&edge.transform(model.datums())?, tol)?;
        let reversed = edge.orientation() == ogeom_topo::Orientation::Reversed;
        let (t0, t1) = if reversed {
            (range.1, range.0)
        } else {
            (range.0, range.1)
        };
        let sign = if reversed { -1.0 } else { 1.0 };
        let unit = |t: f64| -> OgeomResult<Vector> {
            let d = placed.d1_at(t, tol)? * sign;
            if d.magnitude() <= tol.confusion() {
                ogeom_bail!(Construction, "a spine edge has no direction at {t}");
            }
            Ok(d / d.magnitude())
        };
        let station = match &placed {
            Curve::Line(_) => Station {
                from: placed.point_at(t0, tol)?,
                to: placed.point_at(t1, tol)?,
                tangent_in: unit(t0)?,
                tangent_out: unit(t1)?,
                turn: None,
            },
            Curve::Circle(c) => {
                let circle = c.circle();
                let swept = (range.1 - range.0).abs();
                let axis = ogeom_math::Axis {
                    location: circle.centre(),
                    direction: if reversed {
                        -circle.frame().z()
                    } else {
                        circle.frame().z()
                    },
                };
                Station {
                    from: placed.point_at(t0, tol)?,
                    to: placed.point_at(t1, tol)?,
                    tangent_in: unit(t0)?,
                    tangent_out: unit(t1)?,
                    turn: Some((axis, swept)),
                }
            }
            other => ogeom_bail!(
                Construction,
                "a spine runs on straight and circular edges; a {:?} sweeps a \
                 surface this construction does not have",
                other.kind()
            ),
        };
        out.push(station);
    }
    Ok(out)
}

/// The spine's own normal, and the check that it has one.
///
/// Taken from the first turn the spine makes (a corner or an arc) because
/// that is exact, and then measured against every station: a spine that
/// leaves its own plane has no square profile to carry, and says so here
/// rather than by producing a shape nobody asked for.
fn spine_normal(stations: &[Station], tol: Tolerances) -> OgeomResult<Direction> {
    let mut best: Option<(f64, Vector)> = None;
    let mut consider = |a: Vector, b: Vector| {
        let cross = a.cross(b);
        let magnitude = cross.magnitude();
        if magnitude > best.map_or(tol.angular(), |(m, _)| m) {
            best = Some((magnitude, cross / magnitude));
        }
    };
    for (index, station) in stations.iter().enumerate() {
        consider(station.tangent_in, station.tangent_out);
        if index + 1 < stations.len() {
            consider(station.tangent_out, stations[index + 1].tangent_in);
        }
    }
    if stations.len() > 1 {
        consider(
            stations[stations.len() - 1].tangent_out,
            stations[0].tangent_in,
        );
    }
    let Some((_, normal)) = best else {
        ogeom_bail!(
            Construction,
            "a spine that never turns has no plane of its own; give the \
             profile's own orientation a spine with at least one corner or arc"
        );
    };
    for station in stations {
        for tangent in [station.tangent_in, station.tangent_out] {
            if tangent.dot(normal).abs() > tol.angular() {
                ogeom_bail!(
                    Construction,
                    "the spine leaves its own plane; an evolved sweep runs a \
                     planar spine"
                );
            }
        }
        if let Some((axis, _)) = &station.turn
            && axis.direction.vector().cross(normal).magnitude() > tol.angular()
        {
            ogeom_bail!(
                Construction,
                "a spine arc turns about an axis off the spine's own normal"
            );
        }
    }
    Direction::new(normal, tol)
}

/// The profile's plane: a point on it and its normal.
fn profile_plane(model: &Model, profile: &Shape, tol: Tolerances) -> OgeomResult<(Point, Vector)> {
    let mut points: Vec<Point> = Vec::new();
    for edge in explore(model, profile, Filter::OfType(ShapeType::Edge))? {
        let Some(data) = model.node(&edge).and_then(|n| n.data().as_edge()) else {
            continue;
        };
        let Some(EdgeRepr::Curve3d { curve, range, .. }) = data.curve3d() else {
            continue;
        };
        let Some(geometry) = model.geometry().curve(*curve) else {
            ogeom_bail!(Dangling, "curve is not in this model");
        };
        let placed = geometry.transformed(&edge.transform(model.datums())?, tol)?;
        for i in 0..=8 {
            let t = range.0 + (range.1 - range.0) * f64::from(i) / 8.0;
            points.push(placed.point_at(t, tol)?);
        }
    }
    if points.len() < 3 {
        ogeom_bail!(Construction, "a profile needs an extent to sweep");
    }
    let origin = points[0];
    // The widest cross product among the sampled offsets: the plane's normal,
    // taken where it is best conditioned rather than from the first three
    // points that happen to be there.
    let mut best: Option<(f64, Vector)> = None;
    for (i, a) in points.iter().enumerate() {
        for b in points.iter().skip(i + 1) {
            let cross = (*a - origin).cross(*b - origin);
            let magnitude = cross.magnitude();
            if magnitude > best.map_or(tol.confusion(), |(m, _)| m) {
                best = Some((magnitude, cross / magnitude));
            }
        }
    }
    let Some((_, normal)) = best else {
        ogeom_bail!(Construction, "a profile with no area has no plane");
    };
    for p in &points {
        if (*p - origin).dot(normal).abs() > tol.confusion() * 1e2 {
            ogeom_bail!(Construction, "the profile is not planar");
        }
    }
    Ok((origin, normal))
}

/// The rigid motion that carries the profile from the spine's start to a
/// station: a turn about the spine's normal, then a translation.
fn station_transform(
    reference: (Point, Vector),
    at: Point,
    tangent: Vector,
    normal: Direction,
    tol: Tolerances,
) -> OgeomResult<Transform> {
    let (origin, from) = reference;
    let n = normal.vector();
    let angle = from.cross(tangent).dot(n).atan2(from.dot(tangent));
    let turn = if angle.abs() <= tol.angular() {
        Transform::IDENTITY
    } else {
        Transform::rotation(
            ogeom_math::Axis {
                location: origin,
                direction: normal,
            },
            angle,
        )
    };
    Ok(Transform::translation(at - origin) * turn)
}

/// The wedge a corner adds: the profile turned about the corner, through
/// exactly the angle the spine turns there.
///
/// `None` where the spine does not turn; two edges meeting smoothly leave no
/// wedge to fill.
#[allow(clippy::too_many_arguments)]
fn corner_piece(
    model: &mut Model,
    profile: &Shape,
    reference: (Point, Vector),
    corner: Point,
    incoming: Vector,
    outgoing: Vector,
    normal: Direction,
    tol: Tolerances,
) -> OgeomResult<Option<Shape>> {
    let n = normal.vector();
    let angle = incoming
        .cross(outgoing)
        .dot(n)
        .atan2(incoming.dot(outgoing));
    if angle.abs() <= tol.angular() {
        return Ok(None);
    }
    let placed = ogeom_algo::transformed(
        model,
        profile,
        station_transform(reference, corner, incoming, normal, tol)?,
    )?
    .shape;
    let axis = ogeom_math::Axis {
        location: corner,
        direction: if angle > 0.0 { normal } else { -normal },
    };
    let turned = ogeom_algo::make_revolution(model, &placed, axis, angle.abs(), tol);
    match turned {
        Ok(built) => Ok(Some(built.shape)),
        Err(_) => ogeom_bail!(
            NotDone,
            "the profile straddles the spine at a corner, so turning it about \
             that corner sweeps it through itself; there is no revolution for \
             that wedge"
        ),
    }
}

/// Pieces laid end to end, fused: pairs of neighbours first, then pairs of
/// those, so each fuse joins two runs of about equal size where one at a
/// time fused every piece onto everything before it, redoing the whole run
/// so far each time. Neighbours stay neighbours, so every fuse still meets
/// its partner across the shared section. `None` for no pieces.
fn fuse_in_order(
    model: &mut Model,
    pieces: Vec<Shape>,
    tol: Tolerances,
) -> OgeomResult<Option<Shape>> {
    let mut runs = pieces;
    while runs.len() > 1 {
        let mut next = Vec::with_capacity(runs.len().div_ceil(2));
        let mut pending = runs.into_iter();
        while let Some(first) = pending.next() {
            match pending.next() {
                Some(second) => next.push(ogeom_bool::fuse(model, &first, &second, tol)?.shape),
                None => next.push(first),
            }
        }
        runs = next;
    }
    Ok(runs.pop())
}