Skip to main content

brep_kernel/construction/sweep_topology/
sweep.rs

1use super::*;
2
3/// Sweep a CLOSED PLANAR profile loop along a path curve (§3.3/§5.7),
4/// TRANSPLANTING the profile onto the path — see [`SectionPlacement`] for the
5/// other reading, and for why the choice changes the solid.
6///
7/// Reuses the loft builder rather than a bespoke swept surface: the path is
8/// sampled at `STATIONS` uniform stations, a rotation-minimizing frame is
9/// propagated along it (double-reflection RMF, Wang et al. 2008 — the frame
10/// does NOT spin at inflections the way raw Frenet does), a rigidly
11/// transformed copy of the profile is placed at each station (identical
12/// degree/knots/weights, only the control points moved — which guarantees
13/// loft's per-section compatibility), and the stations are lofted through
14/// to produce the tube plus planar end caps.
15///
16/// Guards return a clear `Err` on an open/non-planar profile, a degenerate
17/// path tangent, or a loft failure. A path that bends TIGHTER than the section
18/// reaches is refused by name before the loft (`ρκ ≥ 1`, §4c of
19/// [`sweep_sections_through_samples`]); it used to be out of scope and the
20/// caller's problem, and what the caller got instead was a wall that passed
21/// through itself.
22pub fn sweep_profile_along_path(
23    profile: &[NurbsCurve],
24    path: &NurbsCurve,
25    name: Option<&str>,
26) -> Result<BrepSolid, String> {
27    // 32 stations is the original fixed sampling — golden parity pins the
28    // emitted geometry to it — and it stays the count for every path whose
29    // tangent turns through at most a full turn. A path that turns further
30    // (a coil) gets the helix builder's density; see [`path_turning_stations`].
31    sweep_profile_along_path_stations(
32        profile,
33        path,
34        name,
35        path_turning_stations(path),
36        0.0,
37        None,
38        SectionPlacement::Transplant,
39    )
40}
41
42/// Twisted path sweep (§5.7): identical to [`sweep_profile_along_path`], but
43/// the profile additionally ROTATES about the path tangent, linearly in ARC
44/// LENGTH, from 0 at the sweep start to `twist_angle` radians (right-handed
45/// about the tangent) at the end.  The arc-length fraction comes from the
46/// sampled station polyline, not the raw path parameter, so a non-uniformly
47/// parameterized path still twists uniformly in space.
48///
49/// STATION LAW: 16 stations per quarter turn of twist, floored at the path
50/// sweep's 32 and capped at 1024 (the loft's dense interpolation solve is
51/// O(stations³) per control column — the same cap the helix uses).  At
52/// 16/quarter-turn the inter-station twist step is Δφ ≈ 5.6°, so the cubic
53/// v-interpolation error on a profile point circling at radius r is
54/// ≈ r·Δφ⁴/384 ≈ 2.4·10⁻¹⁰·r — far below any geometric tolerance.  The cap
55/// holds that density up to 16 full turns (|twist| = 32π = 1024/16 quarter
56/// turns); a larger twist would silently alias under the cap, so it is
57/// REJECTED with an honest error instead.  Everything else the path sweep
58/// documents (profile validity, the tight-bend refusal) applies unchanged, and
59/// the tight-bend guard measures the TWISTED section's reach because it reads
60/// the axes the sections were actually placed on; twisting about the profile's
61/// own centroid adds no new radial extent of its own.
62pub fn sweep_profile_twisted(
63    profile: &[NurbsCurve],
64    path: &NurbsCurve,
65    twist_angle: f64,
66    name: Option<&str>,
67) -> Result<BrepSolid, String> {
68    use std::f64::consts::{FRAC_PI_2, TAU};
69
70    if !twist_angle.is_finite() {
71        return Err("sweep_profile_twisted: twist angle must be finite".into());
72    }
73    // 16 turns is where the 1024-station cap meets 16 stations/quarter-turn;
74    // beyond it the cap would degrade the twist sampling density silently.
75    const MAX_TURNS: f64 = 16.0;
76    if twist_angle.abs() > MAX_TURNS * TAU {
77        return Err(format!(
78            "sweep_profile_twisted: twist of {:.3} turns exceeds the {MAX_TURNS}-turn \
79             limit the 1024-station cap can resolve at 16 stations per quarter turn; \
80             split the sweep or reduce the twist",
81            twist_angle.abs() / TAU
82        ));
83    }
84    let quarter_turns = (twist_angle.abs() / FRAC_PI_2).ceil() as usize;
85    let stations = (quarter_turns * 16).clamp(32, 1024).max(path_turning_stations(path));
86    sweep_profile_along_path_stations(
87        profile,
88        path,
89        name,
90        stations,
91        twist_angle,
92        None,
93        SectionPlacement::Transplant,
94    )
95    .map_err(|error| format!("sweep_profile_twisted: {error}"))
96}
97
98/// The station budget of a single-curve path sweep, from how far the path's
99/// TANGENT turns: 8 stations per quarter turn of tangent — the 32 per turn
100/// [`sweep_profile_helix`] lofts a coil with — floored at the original 32 and
101/// capped at the loft's 1024.
102///
103/// A fixed 32 was a budget for a path that bends once. Along a coil it is
104/// FOUR stations a turn (an 8-turn helix edge, 2026-09-25), and a loft through
105/// sections a quarter turn apart does not follow the path between them: the
106/// wall cut its own corner and folded through itself, and the Tube and Path
107/// Sweep features refused the edge. The floor keeps every path that turns
108/// through at most a full turn on exactly the geometry it always had.
109///
110/// The turning is [`segment_turning`]'s — the chain sampler's own reading —
111/// independent of the station count it decides. A path it cannot read is
112/// reported by the sampler itself, in its own words, after the profile has been
113/// checked; here it only leaves the budget on the floor.
114fn path_turning_stations(path: &NurbsCurve) -> usize {
115    use std::f64::consts::FRAC_PI_2;
116    let turning = path
117        .domain()
118        .and_then(|domain| segment_turning(path, domain))
119        .unwrap_or(0.0);
120    // The slack keeps a path of EXACTLY a full turn (a circle) on the floor
121    // rather than one quarter over it from sampling round-off.
122    let quarter_turns = (turning / FRAC_PI_2 - 1e-3).ceil().max(0.0) as usize;
123    (quarter_turns * 8).clamp(32, 1024)
124}
125
126/// HOW a section is placed at each station — the two things "sweep a profile
127/// along a path" can mean, and the reason they are not the same shape.
128///
129/// `Transplant` MOVES the profile onto the path: the anchor origin lands on the
130/// station point and the profile plane becomes the station's NORMAL plane. The
131/// profile's own position and its angle to the path are both discarded — every
132/// sweep comes out centred on the path and square to it. That is the classic
133/// swept-surface reading, and it is what `SWP` builds.
134///
135/// `Rigid` CARRIES the profile: the section at station `k` is the profile moved
136/// by the same rigid motion the PATH undergoes between its start and station
137/// `k`,
138///
139/// ```text
140///     section_k(x) = P_k + R_k · (x − P_0)
141/// ```
142///
143/// where `R_k` is the rotation taking the start frame `(T_0, r_0, s_0)` to the
144/// station's `(T_k, r_k, s_k)`. The profile keeps its drawn POSITION and its
145/// drawn ANGLE to the path; what changes along the sweep is only what the path
146/// itself does. Two consequences are the whole point of the mode:
147///   - a STRAIGHT path makes every `R_k` the identity, so the section merely
148///     translates and the result is the oblique prism `extrude_profile_brep`
149///     builds from the same profile and vector — `pathAlign` and `translate`
150///     agree on a straight path rather than merely having the same volume;
151///   - a circular ARC makes `R_k` the rotation about the ARC'S OWN CENTRE AXIS
152///     (for a planar path the RMF's transport rotation IS that rotation), so
153///     the profile is carried round the pivot and the end cap arrives at the
154///     angle the start cap had to the tangent. A full-circle path degenerates
155///     to exactly a revolve of the profile about that axis;
156///   - a HELIX makes `R_k` the SCREW about the helix's own axis — but only when
157///     the path says it is one ([`SweepPath::screw_axes`]), because the frame
158///     that realizes it is not the RMF. A helix is carried into itself by the
159///     screw, which spins the tangent frame at the helix's torsion; the RMF has
160///     no spin about the tangent by construction, so it rolls away from the
161///     screw by `∫τ ds` (169.75° over three turns of a radius-5, pitch-5 coil).
162///     A helix segment's stations are therefore framed off its axis
163///     (`screw_reference`), and every other segment keeps the RMF — which on a
164///     chain of lines and arcs in ANY planes is still exactly each arc's turn
165///     about its own centre, composed. The two hand over at a joint without a
166///     roll.
167///
168/// Why not one frame for every curve: the Frenet frame is the screw on a helix
169/// but spins where a spatial curve nearly straightens (a planar S-cubic lifted
170/// 1e-3 out of its plane carries a 180° roll through the inflection, which 32
171/// stations read as 0.4°), and one axis fitted to the whole run is the screw on
172/// a helix but rolls a two-arc spatial chain by up to 256° against its arcs.
173///
174/// `R_k = F_k · F_0ᵀ` does NOT depend on which perpendicular `r_0` the RMF was
175/// seeded with: re-seeding rolls `F_0` and every `F_k` by the same angle about
176/// the tangent, and the two rolls cancel in the product. That matters because
177/// `Transplant` maps `pu → r_k` directly, so under IT the profile's roll about
178/// the path is whatever `tangents[0].perpendicular()` happened to return — an
179/// arbitrary orientation no caller can predict. `Rigid` has no such freedom.
180#[derive(Debug, Clone, Copy, PartialEq, Eq)]
181pub enum SectionPlacement {
182    /// Move the profile onto the path, square to it (the classic sweep).
183    Transplant,
184    /// Carry the profile by the path's own rigid motion (see above).
185    Rigid,
186}
187
188/// The placement anchor a path sweep transplants its profile with: the plane
189/// frame `(origin = boundary-sample centroid, normal, pu, pv)` the station
190/// loop maps profile points through (`local = p − origin` → `station + ri·(local
191/// ·pu) + si·(local·pv)`). Extracted as data so a HOLE loop can sweep with its
192/// OUTER loop's anchor — sweeping each loop with its OWN centroid would
193/// re-center every loop onto the path and lose the hole's in-plane offset.
194#[derive(Debug, Clone, Copy)]
195pub struct ProfileAnchor {
196    pub origin: Vec3,
197    pub normal: Vec3,
198    pub pu: Vec3,
199    pub pv: Vec3,
200}
201
202/// [`sweep_profile_along_path`] with an explicit placement anchor (see
203/// [`ProfileAnchor`]) — the hole-loop cutter path: the swept loop is validated
204/// as usual but PLACED in its outer loop's frame.
205pub fn sweep_profile_along_path_anchored(
206    profile: &[NurbsCurve],
207    path: &NurbsCurve,
208    name: Option<&str>,
209    anchor: ProfileAnchor,
210) -> Result<BrepSolid, String> {
211    sweep_profile_along_path_stations(
212        profile,
213        path,
214        name,
215        path_turning_stations(path),
216        0.0,
217        Some(anchor),
218        SectionPlacement::Transplant,
219    )
220}
221
222/// [`sweep_profile_twisted`] with an explicit placement anchor: the hole loop
223/// twists about the SAME path axis as its outer loop (shared anchor), so the
224/// cutter stays registered with the outer wall through the whole twist.
225pub fn sweep_profile_twisted_anchored(
226    profile: &[NurbsCurve],
227    path: &NurbsCurve,
228    twist_angle: f64,
229    name: Option<&str>,
230    anchor: ProfileAnchor,
231) -> Result<BrepSolid, String> {
232    use std::f64::consts::{FRAC_PI_2, TAU};
233
234    if !twist_angle.is_finite() {
235        return Err("sweep_profile_twisted: twist angle must be finite".into());
236    }
237    const MAX_TURNS: f64 = 16.0;
238    if twist_angle.abs() > MAX_TURNS * TAU {
239        return Err(format!(
240            "sweep_profile_twisted: twist of {:.3} turns exceeds the {MAX_TURNS}-turn \
241             limit the 1024-station cap can resolve at 16 stations per quarter turn; \
242             split the sweep or reduce the twist",
243            twist_angle.abs() / TAU
244        ));
245    }
246    let quarter_turns = (twist_angle.abs() / FRAC_PI_2).ceil() as usize;
247    let stations = (quarter_turns * 16).clamp(32, 1024).max(path_turning_stations(path));
248    sweep_profile_along_path_stations(
249        profile,
250        path,
251        name,
252        stations,
253        twist_angle,
254        Some(anchor),
255        SectionPlacement::Transplant,
256    )
257    .map_err(|error| format!("sweep_profile_twisted: {error}"))
258}
259
260/// Sweep a CLOSED PLANAR profile along a CHAIN of path curves joined
261/// head-to-tail — the multi-segment entry point, with the twist and anchor
262/// options the single-curve family exposes as separate functions folded into
263/// one signature (the feature layer picks all three per loop, so splitting them
264/// four ways here would only push the same match into the caller).
265///
266/// `path` carries the chain AND its classification ([`SweepPath`]): the segment
267/// names for the joint refusal (a short or empty name list degrades to
268/// `<unnamed>` rather than failing), whether the run is CLOSED, the plane it lies
269/// in if it lies in one, and the continuity of every joint. The refusals below
270/// READ that classification rather than each re-deriving it — the corner gate no
271/// longer measures its own tangent break while sampling, and closure is a
272/// property of the path rather than something the sampler notices about its first
273/// and last station.
274///
275/// `corner_advice` is appended to the cornered-joint refusal and is the CALLER's
276/// to write, because the builder has two callers with different answers: path
277/// sweep sends the user to SW (a different feature), while SW's own `pathAlign`
278/// sends them to its `translate` mode. A builder that named one feature would be
279/// telling half its users to switch to the feature they are already in.
280///
281/// A ONE-segment OPEN chain delegates to the single-curve builders unchanged, so
282/// the geometry every existing path sweep emits is untouched by this entry point
283/// existing — the multi-segment sampler is reached only by a path that actually
284/// has a joint. A one-segment CLOSED chain (a resident circular edge) does NOT
285/// delegate: the single-curve samplers cap both ends, which is exactly what a
286/// ring must not do.
287///
288/// A CORNER (a joint past [`MAX_JOINT_TANGENT_BREAK`]) selects the MITRE lane
289/// when every segment is straight — see `sweep_topology/miter.rs` for the
290/// construction and [`CornerPolicy`] for why this entry mitres while
291/// [`sweep_profile_along_chain_with_stations`] refuses. A CLOSED cornered
292/// polyline is a picture FRAME: its closing joint is a mitre like the rest and
293/// the result has no caps. A corner at a CURVED segment is refused by name.
294pub fn sweep_profile_along_chain(
295    profile: &[NurbsCurve],
296    path: &SweepPath,
297    twist_angle: f64,
298    name: Option<&str>,
299    anchor: Option<ProfileAnchor>,
300    placement_mode: SectionPlacement,
301    corner_advice: &str,
302) -> Result<BrepSolid, String> {
303    use std::f64::consts::{FRAC_PI_2, TAU};
304
305    if path.is_empty() {
306        return Err("sweepSolid: path chain is empty".into());
307    }
308    // A bend TIGHTER than the section is the envelope lane's, not the loft's, and
309    // it is asked BEFORE the single-curve shortcut below — that shortcut lands in
310    // the same skinning builder, so taking it first would hide the one path shape
311    // this lane exists for (one circular arc) behind it.
312    let envelope_lane = super::envelope::recognize_swept_envelope(profile, path, placement_mode)?
313        .is_ok();
314
315    // Single segment under `Transplant`: the existing builders, byte for byte.
316    // `Rigid` has no single-curve wrapper to delegate to and falls through to the
317    // chain sampler, which samples a 1-element chain over exactly the domain the
318    // single-curve sampler would (no joints, the whole budget on one segment).
319    // A CLOSED single segment falls through too — the single-curve samplers end
320    // in a capped loft, and a ring has no caps.
321    if let ([single], SectionPlacement::Transplant, false, false) =
322        (path.curves.as_slice(), placement_mode, path.closed, envelope_lane)
323    {
324        return match (anchor, twist_angle == 0.0) {
325            (None, true) => sweep_profile_along_path(profile, single, name),
326            (Some(anchor), true) => {
327                sweep_profile_along_path_anchored(profile, single, name, anchor)
328            }
329            (None, false) => sweep_profile_twisted(profile, single, twist_angle, name),
330            (Some(anchor), false) => {
331                sweep_profile_twisted_anchored(profile, single, twist_angle, name, anchor)
332            }
333        };
334    }
335
336    // Everything else — any multi-segment chain, and a single segment under
337    // `Rigid` — takes the same station law the single-curve twisted builder
338    // uses, so density scales with twist identically on a chained path.
339    if !twist_angle.is_finite() {
340        return Err("sweep_profile_along_chain: twist angle must be finite".into());
341    }
342    const MAX_TURNS: f64 = 16.0;
343    if twist_angle.abs() > MAX_TURNS * TAU {
344        return Err(format!(
345            "sweep_profile_along_chain: twist of {:.3} turns exceeds the {MAX_TURNS}-turn \
346             limit the 1024-station cap can resolve at 16 stations per quarter turn; \
347             split the sweep or reduce the twist",
348            twist_angle.abs() / TAU
349        ));
350    }
351    // A CLOSED SPATIAL path lays its counter-twist down on top of the requested
352    // twist, which adds up to half a turn, so the density is set for the most the
353    // lap can roll.
354    let stations = if twist_angle == 0.0 {
355        32
356    } else {
357        let counter = if path.closed && path.planar.is_none() {
358            std::f64::consts::PI
359        } else {
360            0.0
361        };
362        let quarter_turns = ((twist_angle.abs() + counter) / FRAC_PI_2).ceil() as usize;
363        (quarter_turns * 16).clamp(32, 1024)
364    };
365    let _ = name;
366    sweep_profile_along_chain_stations(
367        profile,
368        path,
369        stations,
370        twist_angle,
371        anchor,
372        placement_mode,
373        corner_advice,
374        CornerPolicy::Mitre,
375    )
376}
377
378/// [`sweep_profile_along_chain`] with an explicit STATION budget and no twist
379/// — the entry a caller with a long, smooth chain uses when the default 32
380/// stations would under-sample it (a wire-harness bundle along a many-anchor
381/// spline: three pieces per span, so a handful of spans already starves each
382/// piece of stations). The budget is spread over the chain by arc length
383/// exactly as the twisted builder's law spreads its own; the 1024-station cap
384/// and the tangent-continuity refusal apply unchanged.
385pub fn sweep_profile_along_chain_with_stations(
386    profile: &[NurbsCurve],
387    path: &SweepPath,
388    stations: usize,
389    corner_advice: &str,
390) -> Result<BrepSolid, String> {
391    sweep_profile_along_chain_stations(
392        profile,
393        path,
394        stations.clamp(2, 1024),
395        0.0,
396        None,
397        SectionPlacement::Transplant,
398        corner_advice,
399        CornerPolicy::Refuse,
400    )
401}
402
403/// Validate a closed planar profile loop and derive its placement anchor —
404/// the path sweep's §1 block, extracted bit-identically: 16 samples per curve,
405/// closure at `tolerance`, Newell normal, boundary-sample-mean origin,
406/// planarity at `tolerance * 100`, `pu = np.perpendicular()`, `pv = np × pu`.
407pub fn profile_anchor(profile: &[NurbsCurve]) -> Result<ProfileAnchor, String> {
408    let tolerance = 1e-6;
409    if profile.len() < 2 {
410        return Err("sweepSolid: profile needs at least 2 curves forming a closed loop".into());
411    }
412    let mut samples = Vec::new();
413    for (index, curve) in profile.iter().enumerate() {
414        let [start, end] = curve.domain()?;
415        let next = &profile[(index + 1) % profile.len()];
416        let next_start = next.domain()?[0];
417        if curve
418            .evaluate(end)?
419            .sub(next.evaluate(next_start)?)
420            .length()
421            > tolerance
422        {
423            return Err(format!(
424                "sweepSolid: profile is not closed at curve {index}"
425            ));
426        }
427        for sample in 0..16 {
428            samples.push(curve.evaluate(start + (end - start) * sample as f64 / 16.0)?);
429        }
430    }
431    let normal = crate::polygon::newell_normal(&samples);
432    let centroid = samples.iter().fold(Vec3::default(), |sum, &point| sum.add(point));
433    let np = normal
434        .normalized()
435        .map_err(|_| "sweepSolid: profile is degenerate (zero enclosed area)".to_string())?;
436    let origin = centroid.scale(1.0 / samples.len() as f64);
437    if samples
438        .iter()
439        .any(|point| point.sub(origin).dot(np).abs() > tolerance * 100.0)
440    {
441        return Err("sweepSolid: profile is not planar".into());
442    }
443    let pu = np.perpendicular()?;
444    let pv = np.cross(pu).normalized()?;
445    Ok(ProfileAnchor {
446        origin,
447        normal: np,
448        pu,
449        pv,
450    })
451}
452
453/// Core of the path sweep with an explicit station count.  Every consumer of
454/// the station count (path sampling, RMF propagation, section placement)
455/// derives from the one `stations` argument so the density scales as a unit.
456/// `twist_angle` (radians; 0 for the untwisted variants) rotates the placed
457/// profile about the path tangent linearly in sampled arc length — the
458/// `twist_angle == 0.0` fast path leaves the RMF axes bit-identical, so the
459/// untwisted callers keep golden parity.
460fn sweep_profile_along_path_stations(
461    profile: &[NurbsCurve],
462    path: &NurbsCurve,
463    name: Option<&str>,
464    stations: usize,
465    twist_angle: f64,
466    anchor: Option<ProfileAnchor>,
467    placement_mode: SectionPlacement,
468) -> Result<BrepSolid, String> {
469    // Loft carries no face names; the app stamps them onto the emitted face
470    // order.  Accept `name` for ABI symmetry with the other builders.
471    let _ = name;
472    let tolerance = 1e-6;
473    if stations < 2 {
474        return Err("sweepSolid: need at least 2 stations".into());
475    }
476
477    // --- 1. Validate the profile: closed + planar; derive (origin, np, pu, pv).
478    //        A caller-supplied anchor OVERRIDES the placement frame (the swept
479    //        loop is still validated against its own plane), so a hole loop
480    //        rides the path in its outer loop's frame instead of re-centering.
481    //        Resolved HERE, before the path is touched, so a bad profile is
482    //        still reported ahead of a bad path exactly as it always was.
483    let placement = resolve_placement(profile, anchor, placement_mode)?;
484
485    // --- 2. Sample the path; require a non-degenerate tangent at every station.
486    let [t0, t1] = path.domain()?;
487    if (t1 - t0).abs() <= tolerance {
488        return Err("sweepSolid: path domain is degenerate".into());
489    }
490    let mut points = Vec::with_capacity(stations);
491    let mut tangents = Vec::with_capacity(stations);
492    let mut curvatures = Vec::with_capacity(stations);
493    for index in 0..stations {
494        let t = t0 + (t1 - t0) * index as f64 / (stations - 1) as f64;
495        let derivatives = path.derivatives_small(t, 2)?;
496        let tangent = derivatives[1]
497            .normalized()
498            .map_err(|_| format!("sweepSolid: path tangent is degenerate at station {index}"))?;
499        points.push(derivatives[0]);
500        tangents.push(tangent);
501        curvatures.push(curvature_vector(derivatives[1], derivatives[2]));
502    }
503
504    sweep_sections_through_samples(
505        profile,
506        &points,
507        &tangents,
508        &curvatures,
509        twist_angle,
510        placement,
511        placement_mode,
512        RingFrame::Open,
513        None,
514        None,
515    )
516}
517
518/// What a CHAIN builder does with a CORNER (a `G0` joint).
519///
520/// Two of this module's entry points sweep chains, and they answer the question
521/// differently — not because the geometry differs but because their callers'
522/// corner POLICIES do, and a builder that guessed would be overriding one of
523/// them:
524///   - [`sweep_profile_along_chain`] MITRES a polyline's corners, a closed
525///     polyline's closing joint included (`sweep_topology/miter.rs`). Its callers
526///     are the two sweep FEATURES, and a mitred corner is the shape they are
527///     asked for.
528///   - [`sweep_profile_along_chain_with_stations`] REFUSES every corner. Its
529///     callers build their own chains and have already decided what a corner
530///     means: the tube ROUNDS one into a tangent bend or joins it with a ball,
531///     and the wire harness routes tangent-continuous splines. A corner
532///     reaching either of them means their own construction did not come out,
533///     and mitring it instead would silently ship a shape neither feature
534///     offers.
535#[derive(Debug, Clone, Copy, PartialEq, Eq)]
536enum CornerPolicy {
537    /// A polyline's corners are mitred, open or closed; anything else is refused.
538    Mitre,
539    /// Every corner is refused by name.
540    Refuse,
541}
542
543/// The SECTION a mitred sweep carries, placed in world, with its plane normal —
544/// the two things `sweep_topology/miter.rs` needs, and all it needs.
545///
546/// The placement modes are resolved HERE rather than in the mitre builder, on
547/// the same frame the skinning lane's station 0 uses (`û₀.perpendicular()`), so
548/// the two lanes place a section identically and a path that crosses the
549/// continuity gate does not also jump in space:
550///   - `Transplant` maps the profile's plane frame onto the path's first station
551///     (anchor origin onto the path, profile square to it), so the section's
552///     normal is the first tangent;
553///   - `Rigid` takes the profile WHERE IT IS, keeping the offset from the path
554///     and the angle to it that it was drawn with, so the section's normal is
555///     the profile's own.
556fn place_section(
557    profile: &[NurbsCurve],
558    placement: ProfileAnchor,
559    placement_mode: SectionPlacement,
560    start: Vec3,
561    tangent: Vec3,
562) -> Result<(Vec<NurbsCurve>, Vec3), String> {
563    if placement_mode == SectionPlacement::Rigid {
564        return Ok((profile.to_vec(), placement.normal));
565    }
566    let ProfileAnchor { origin, pu, pv, .. } = placement;
567    let r0 = tangent.perpendicular()?;
568    let s0 = tangent.cross(r0).normalized()?;
569    let mut placed = Vec::with_capacity(profile.len());
570    for curve in profile {
571        let control_points = curve
572            .control_points
573            .iter()
574            .map(|control| {
575                let point = control.point()?;
576                let local = point.sub(origin);
577                let world = start
578                    .add(r0.scale(local.dot(pu)))
579                    .add(s0.scale(local.dot(pv)));
580                Ok(Vec4::from_point(world, control.w))
581            })
582            .collect::<Result<Vec<_>, String>>()?;
583        placed.push(NurbsCurve::new(
584            curve.degree,
585            curve.knots.clone(),
586            control_points,
587        )?);
588    }
589    Ok((placed, tangent))
590}
591
592/// The placement frame a sweep transplants its profile with: the caller's
593/// `anchor` when it supplied one (a hole loop riding its outer loop's frame),
594/// otherwise the profile's own. Either way the profile is VALIDATED — closed and
595/// planar — because `profile_anchor` is what does that validating, and skipping
596/// it for an anchored loop would let an open hole loop through.
597///
598/// Called by each SAMPLER before it touches the path, so the "bad profile" error
599/// still precedes the "bad path" one, as it did when the two lived in one
600/// function.
601fn resolve_placement(
602    profile: &[NurbsCurve],
603    anchor: Option<ProfileAnchor>,
604    placement_mode: SectionPlacement,
605) -> Result<ProfileAnchor, String> {
606    let computed = profile_anchor(profile)?;
607    Ok(match placement_mode {
608        // RIGID never transplants, so there is no frame to borrow: every loop
609        // keeps its own position already, which is the whole reason the anchor
610        // exists under `Transplant`. Taking a caller's anchor here would silently
611        // swap in another loop's NORMAL for the fold-back guard below.
612        SectionPlacement::Rigid => computed,
613        SectionPlacement::Transplant => anchor.unwrap_or(computed),
614    })
615}
616
617
618/// Sweep a CLOSED PLANAR profile along a CHAIN of path curves, joined
619/// head-to-tail, as one continuous tube — or, for a CLOSED path, as a capless
620/// ring.
621///
622/// The chain is sampled as a single trajectory and handed to the shared core,
623/// so ONE rotation-minimizing frame is propagated across the whole path: a joint
624/// is just another pair of adjacent stations to the RMF, which is what makes the
625/// profile arrive at the far end in the orientation the near end implies rather
626/// than snapping at each segment.
627///
628/// STATION BUDGET. `stations` are distributed across the segments in proportion
629/// to estimated ARC LENGTH (chord sums over 16 samples per curve), floored at 2
630/// per segment so a short fillet between two long lines still carries its
631/// curvature. The shared endpoint of two adjacent segments is sampled ONCE, not
632/// twice — a duplicated station would land in the RMF's coincident-step branch
633/// and waste a section on zero advance. A CLOSED path drops one more station for
634/// the same reason: its last station IS its first.
635///
636/// JOINTS decide the LANE. A break within [`MAX_JOINT_TANGENT_BREAK`] is skinned
637/// here; a break past it is a CORNER, which under [`CornerPolicy::Mitre`] hands an
638/// open all-straight run to the mitre builder (`sweep_topology/miter.rs`) and is
639/// otherwise refused by name, with the angle and the two segments. The
640/// classification carried on the path is what says so — this sampler no longer
641/// measures its own.
642#[allow(clippy::too_many_arguments)]
643fn sweep_profile_along_chain_stations(
644    profile: &[NurbsCurve],
645    path: &SweepPath,
646    stations: usize,
647    twist_angle: f64,
648    anchor: Option<ProfileAnchor>,
649    placement_mode: SectionPlacement,
650    corner_advice: &str,
651    corners: CornerPolicy,
652) -> Result<BrepSolid, String> {
653    let chain = path.curves.as_slice();
654    if chain.is_empty() {
655        return Err("sweepSolid: path chain is empty".into());
656    }
657    if stations < 2 {
658        return Err("sweepSolid: need at least 2 stations".into());
659    }
660    // Profile first, matching the single-curve sampler's order of complaint.
661    let placement = resolve_placement(profile, anchor, placement_mode)?;
662
663    // --- The swept ENVELOPE, trimmed at its own self-intersection.
664    //
665    //     Past `ρκ = 1` the sections CROSS and the loft below skins a wall that
666    //     passes through itself. For a CIRCLE carried along one circular planar
667    //     arc the correct solid is a REVOLVE — the section disc truncated at the
668    //     axis, plus the reflection of the part that crossed it, each revolved
669    //     through the path's own turn (`sweep_topology/envelope.rs`) — and this
670    //     is where that lane is taken instead of the loft.
671    //
672    //     Anything else tight is left to the §4c refusal below, which now carries
673    //     WHICH half of the configuration this sweep missed. Nothing under
674    //     `ρκ = 1` reaches the lane at all: the recognition refuses a bend that
675    //     is not tighter than its section, so every shape that builds today keeps
676    //     the builder it has.
677    let envelope_gap =
678        match super::envelope::recognize_swept_envelope(profile, path, placement_mode)? {
679            Ok(plan) => {
680                let mut bodies = super::envelope::build_swept_envelope(&plan)?;
681                if bodies.len() == 1 {
682                    return Ok(bodies.remove(0));
683                }
684                return Err(format!(
685                    "{}: a section that crosses the path's axis (ρ·κ = {:.6}) sweeps TWO \
686                     bodies — the near lobe over the path's own azimuths and the far lobe half \
687                     a turn from them, meeting only along the axis — and this entry returns one \
688                     solid. {} were built; ask `sweep_envelope_bodies` for the pair",
689                    SWEEP_ENVELOPE_REFUSAL,
690                    plan.ratio(),
691                    bodies.len()
692                ));
693            }
694            Err(reason) => reason,
695        };
696
697    // --- A CORNER, read off the path's classification. Every joint was measured
698    //     once when the path was resolved, the CLOSING joint of a closed path
699    //     included — which is how a ring's seam gets the same gate as its
700    //     interior joints instead of being the one joint nobody checked.
701    //
702    //     A corner now selects a LANE rather than ending the build. A POLYLINE —
703    //     every segment straight, at least one corner — is MITRED: the section is
704    //     swept along each segment and the two sweeps are trimmed to the joint's
705    //     bisector plane, sharing one face loop there
706    //     (`sweep_topology/miter.rs`). A CLOSED one is a picture FRAME, whose
707    //     closing joint is mitred like the rest and which has no caps; the frame
708    //     lane's own refusals (a spatial loop, by its holonomy; a self-crossing
709    //     loop; overlapping mitres) live with the construction. What remains
710    //     refused HERE is a corner on a path that is not all straight — there is
711    //     no bisector plane for a curve to be trimmed to, and no lane that mixes
712    //     the two.
713    if let Some(corner) = path.corner() {
714        if corners == CornerPolicy::Mitre && path.cornered_polyline() {
715            // A twist is laid down in arc length over the stations, and a mitre has
716            // none. A CLOSED frame rolls it along the middle of each side instead,
717            // where a roll keeps the shared joint loops (`sweep_topology/miter.rs`),
718            // exactly as it rolls its own counter-twist. An OPEN cornered run does
719            // not have that construction and is refused.
720            if twist_angle != 0.0 && !path.closed {
721                return Err(
722                    "sweepSolid: a CORNERED path cannot carry a twist — its corners are mitred \
723                     from exact pieces that meet on the joint's bisector plane, and a twist laid \
724                     down along the path would roll the two sides of that plane apart, so they \
725                     would no longer share one loop there. Only a CLOSED frame rolls a twist along \
726                     the middle of its sides. Sweep the run without the twist; a \
727                     tangent-continuous (smooth) path is skinned through stations and carries one"
728                        .into(),
729                );
730            }
731            let (section, section_normal, _) =
732                mitre_section(profile, placement, placement_mode, path)?;
733            return super::miter::sweep_section_mitered(
734                &section,
735                section_normal,
736                path,
737                twist_angle,
738            );
739        }
740        // A path with BOTH a corner and a curved segment is in NEITHER lane, and
741        // this is the refusal that says so: the mitre lane needs every segment
742        // straight (one direction per side for the bisector plane to bisect) and
743        // the skinning lane needs every joint tangent-continuous (it interpolates
744        // stations, so it would round the corner off). A rounded rectangle with
745        // ONE corner left sharp is exactly that path — three G0 joints the mitre
746        // lane could carry and two G1 joints at an arc it cannot — and mitring
747        // the sharp corners while skinning the arc is a MIXED construction this
748        // builder does not have. Refused with both halves named rather than built
749        // as whichever lane happens to be asked first.
750        let curved = (0..path.len())
751            .find(|index| {
752                let curve = &path.curves[*index];
753                curve.degree != 1 || curve.control_points.len() != 2
754            })
755            .map(|index| format!(" — segment '{}' is CURVED", path.name(index)))
756            .unwrap_or_default();
757        // The MIXED-LANE sentence belongs only to a path that actually has a
758        // curved segment. Under [`CornerPolicy::Refuse`] an ALL-STRAIGHT cornered
759        // chain reaches this same refusal (the tube's and the wire harness's
760        // corners are refused whatever their segments are), and telling those
761        // callers that a straight cornered run is mitred would name a lane they
762        // deliberately do not have — their own `corner_advice` says what to do
763        // instead.
764        let mixed = if curved.is_empty() {
765            String::new()
766        } else {
767            " A path that MIXES the two — corners at some joints and curved segments elsewhere              — is in neither lane: every segment straight is mitred (open or closed), every              joint tangent-continuous is skinned, and nothing builds one construction for some              joints and the other for the rest."
768                .to_string()
769        };
770        return Err(format!(
771            "the path must be tangent-continuous: segments '{}' and '{}' meet at a {:.1}° \
772             corner{curved}. The profile is skinned between sampled stations along the path, so a \
773             corner would be ROUNDED OFF rather than mitred, and a mitre needs two STRAIGHT \
774             segments: it trims both sweeps to the plane that bisects their directions, and a \
775             curve has no single direction to bisect.{mixed} {corner_advice}",
776            path.name(corner.before),
777            path.name(corner.after),
778            corner.tangent_break.to_degrees()
779        ));
780    }
781
782    // --- A CLOSED path takes a REQUESTED twist only where it closes. After one lap
783    //     the section comes back rolled by the whole twist, so the twist must be a
784    //     rotational symmetry of the section about the path. The ring lane reads
785    //     that off the section and refuses by name a twist that does not close
786    //     ([`closing_twist`]).
787
788    // --- Sample the chain as ONE trajectory (budget, domains, stations).
789    let StationSamples {
790        points,
791        tangents,
792        curvatures,
793        segments,
794    } = chain_stations(path, stations)?;
795
796    // --- CLOSED: the ring needs a frame that comes back to itself. A PLANAR path
797    //     has one for free (see `ring_normal`); a SPATIAL one comes back rolled by
798    //     its own holonomy and is closed by the counter-twist policy
799    //     ([`RingFrame::Spatial`]).
800    let closure = ring_frame(path, &tangents)?;
801
802    sweep_sections_through_samples(
803        profile,
804        &points,
805        &tangents,
806        &curvatures,
807        twist_angle,
808        placement,
809        placement_mode,
810        closure,
811        Some(envelope_gap.as_str()),
812        Some((segments.as_slice(), path)),
813    )
814}
815
816/// How a station run CLOSES — which frame a ring's sections are placed on.
817///
818/// Three states, because "is this a ring" and "is its frame written down in
819/// closed form" are different questions and only a planar ring answers yes to
820/// both:
821///   - `Open`: the transported frame and two end caps;
822///   - `Planar(normal)`: a ring whose reference axis is the path plane's normal
823///     at every station, which returns to itself EXACTLY (see [`ring_normal`]);
824///   - `Spatial`: a ring on the TRANSPORTED frame, which returns rolled by the
825///     path's own holonomy and is closed by the counter-twist policy — see
826///     [`ring_closure_frames`].
827#[derive(Debug, Clone, Copy)]
828enum RingFrame {
829    Open,
830    Planar(Vec3),
831    Spatial,
832}
833
834impl RingFrame {
835    fn is_ring(self) -> bool {
836        !matches!(self, RingFrame::Open)
837    }
838}
839
840/// The ring frame a path's stations close on, read off its classification.
841fn ring_frame(path: &SweepPath, tangents: &[Vec3]) -> Result<RingFrame, String> {
842    Ok(match (path.closed, path.planar) {
843        (false, _) => RingFrame::Open,
844        (true, Some(normal)) => RingFrame::Planar(ring_normal(normal, tangents)?),
845        (true, None) => RingFrame::Spatial,
846    })
847}
848
849/// What a CLOSED sweep did to make its section meet itself after one lap — the
850/// applied twist the frame-closure policy REPORTS.
851///
852/// THE POLICY. Carried once round a closed path, the section frame comes back
853/// ROLLED about the start tangent by the path's HOLONOMY. That angle is not a
854/// defect of the sampling: the frame is parallel transport of the section's
855/// axes in the path's normal plane, which is parallel transport on the unit
856/// sphere along the path's TANGENT INDICATRIX, so by Gauss–Bonnet the roll is
857/// the solid angle that indicatrix encloses, modulo `2π`. A planar path's
858/// indicatrix runs along one great circle and encloses a hemisphere — `2π`,
859/// the identity — which is why a planar ring closes on its own; a spatial one
860/// generally does not. (A curve lying ON a sphere reads zero too: the radial
861/// direction is itself rotation-minimizing along any such curve. Its enclosed
862/// solid angle is not its holonomy; its indicatrix's is.)
863///
864/// The sweep closes the ring with the SMALLEST roll that does it: the
865/// counter-twist `−holonomy`, wrapped to `(−π, π]`.
866///   - SMOOTH lane: laid down LINEARLY IN ARC LENGTH over the stations, the wrap
867///     step from the last station back to station 0 included, so every span —
868///     the seam span too — rolls at the same rate and the closed loft sees no
869///     seam at all. A uniform rate is the distribution with the smallest peak
870///     rate and the smallest `∫(roll rate)² ds`.
871///   - CORNERED lane (`sweep_topology/miter.rs`): shared over the SEGMENTS in
872///     proportion to their length, each share laid down over its segment's
873///     untrimmed middle — a roll AT a mitre would pull the shared loop apart —
874///     and nothing at all when the untwisted lap already closes within the
875///     lane's own weld floor.
876///
877/// A REQUESTED twist rides on top of that counter-twist and is laid down the same
878/// way, so `applied_twist = requested_twist − holonomy` (or `requested_twist` alone
879/// where no counter-twist is needed). It closes only as a rotational symmetry of
880/// the section about the path, `2π·n/k` for a section with `k` of them, and is
881/// refused by name otherwise ([`crate::SWEEP_TWIST_CLOSURE_REFUSAL`]). A
882/// FRACTIONAL turn on a symmetric section lands each wall on a different section
883/// curve than it left: `curve_shift`.
884///
885/// `holonomy` is what the untwisted frame measures; `applied_twist` is what was
886/// laid down. All three angles are radians, right-handed about the path tangent.
887#[derive(Debug, Clone, Copy, PartialEq)]
888pub struct SweepClosure {
889    /// The roll the UNTWISTED frame returns with after one lap, right-handed
890    /// about the start tangent, wrapped to `(−π, π]`. Zero for a planar path.
891    pub holonomy: f64,
892    /// The twist the caller asked for, snapped to the section symmetry it closes
893    /// on (it matched at round-off). Zero for an untwisted sweep.
894    pub requested_twist: f64,
895    /// The twist laid down over the lap: the requested twist plus the
896    /// counter-twist `−holonomy`, which is zero where the lane closes without one
897    /// (a planar path, or a frame whose untwisted lap already closes within the
898    /// weld floor).
899    pub applied_twist: f64,
900    /// How far along the section's curve list one lap carries each curve: the
901    /// wall leaving section curve `j` arrives on curve `(j + curve_shift) mod m`.
902    /// Zero for no twist and for a whole number of turns. The curves are counted
903    /// in the order the builder carries them, which on the MITRE lane is the
904    /// section's order after its winding is made counter-clockwise about the path.
905    pub curve_shift: usize,
906}
907
908/// The CLOSURE a sweep of `profile` along `path` applies — the frame-closure
909/// policy's report, exposed so a caller can read the twist a closed sweep laid
910/// down (see [`SweepClosure`]). `None` for an OPEN path, which has no lap to
911/// close.
912///
913/// Reads the very construction the sweep builds, not a second measurement of
914/// it: a cornered polyline asks the mitre lane's own frame plan (its placed
915/// section, its rings, its weld floor), and any other closed path asks the
916/// smooth lane's own stations and frames at `stations` — the budget the sweep
917/// was given (32 for [`sweep_profile_along_chain`] untwisted, the caller's for
918/// [`sweep_profile_along_chain_with_stations`]). Whatever the builder would
919/// refuse, this refuses with the same words, a requested `twist_angle` that does
920/// not close included.
921pub fn sweep_closure(
922    profile: &[NurbsCurve],
923    path: &SweepPath,
924    placement_mode: SectionPlacement,
925    stations: usize,
926    twist_angle: f64,
927) -> Result<Option<SweepClosure>, String> {
928    if !path.closed {
929        return Ok(None);
930    }
931    let placement = resolve_placement(profile, None, placement_mode)?;
932    if path.cornered_polyline() {
933        let (section, section_normal, _) = mitre_section(profile, placement, placement_mode, path)?;
934        return super::miter::mitred_frame_closure(&section, section_normal, path, twist_angle)
935            .map(Some);
936    }
937    let StationSamples {
938        points,
939        tangents,
940        segments,
941        ..
942    } = chain_stations(path, stations.clamp(2, 1024))?;
943    let screw = Some((segments.as_slice(), path)).filter(|_| placement_mode == SectionPlacement::Rigid);
944    let closure = ring_frame(path, &tangents)?;
945    let (section, _) = place_section(profile, placement, placement_mode, points[0], tangents[0])?;
946    let requested = ring_closing_twist(&section, &points, &tangents, twist_angle)?;
947    Ok(station_frames(&points, &tangents, closure, screw, requested)?.2)
948}
949
950/// The requested twist a RING closes on ([`closing_twist`]), read off the
951/// section as placed at station 0, about the path's first tangent, against the
952/// ring's own position floor: `1e-6` of the lap's chord length, the floor the
953/// mitred frame measures its own closure against.
954fn ring_closing_twist(
955    section: &[NurbsCurve],
956    points: &[Vec3],
957    tangents: &[Vec3],
958    twist_angle: f64,
959) -> Result<ClosingTwist, String> {
960    let lap: f64 = (0..points.len())
961        .map(|index| points[(index + 1) % points.len()].sub(points[index]).length())
962        .sum();
963    closing_twist(section, points[0], tangents[0], twist_angle, 1e-6 * lap)
964}
965
966/// The placed section a MITRED sweep carries and its plane normal, with the
967/// first segment's chord direction — shared by the builder and by
968/// [`sweep_closure`], so the instrument reads the section the frame is built
969/// from.
970fn mitre_section(
971    profile: &[NurbsCurve],
972    placement: ProfileAnchor,
973    placement_mode: SectionPlacement,
974    path: &SweepPath,
975) -> Result<(Vec<NurbsCurve>, Vec3, Vec3), String> {
976    let [t0, t1] = path.curves[0].domain()?;
977    let start = path.curves[0].evaluate(t0)?;
978    // The CHORD direction, which is both the first station's tangent (a
979    // straight segment's derivative is its chord) and the direction the
980    // mitre builder sweeps along, so the placed section is square to the
981    // same vector the construction uses.
982    let tangent = path.curves[0]
983        .evaluate(t1)?
984        .sub(start)
985        .normalized()
986        .map_err(|_| "sweepSolid: path segment 0 has zero length".to_string())?;
987    let (section, section_normal) =
988        place_section(profile, placement, placement_mode, start, tangent)?;
989    Ok((section, section_normal, tangent))
990}
991
992/// The station frames §3 places its sections on, extracted so the tight-bend
993/// INSTRUMENT ([`sweep_bend_profile`]) reads the frame the sweep will use — the
994/// section's reach into a bend depends on how the section is ROLLED about the
995/// path, so a second frame that merely agreed on a circle would be a different
996/// measurement on anything else.
997///
998/// `screw` is the chain's segment-of-station map and the path whose
999/// [`SweepPath::screw_axes`] it indexes, passed only by the `Rigid` placement.
1000/// A station on a HELIX segment is framed off the helix's own axis instead of
1001/// transported (see [`screw_reference`]); every other station takes the RMF step
1002/// from its predecessor, so the two compose across a joint in either direction.
1003///
1004/// The third value is the ring's CLOSURE report: `None` for an open run, zero
1005/// for a planar ring, and the measured holonomy with the counter-twist that
1006/// closes it for a spatial one — whose axes are returned already twisted.
1007///
1008/// `requested` is a RING's requested twist, already known to close
1009/// ([`closing_twist`]). It is laid down in the ring's frames over the whole lap,
1010/// the wrap chord included, together with any counter-twist, so the seam span
1011/// rolls at the rate every other span does. An OPEN run ignores it: its twist is
1012/// laid down by the placement, over the stations it has.
1013fn station_frames(
1014    points: &[Vec3],
1015    tangents: &[Vec3],
1016    closure: RingFrame,
1017    screw: Option<(&[usize], &SweepPath)>,
1018    requested: ClosingTwist,
1019) -> Result<(Vec<Vec3>, Vec<Vec3>, Option<SweepClosure>), String> {
1020    match closure {
1021        RingFrame::Open => {
1022            let (r_axes, s_axes) = transported_frames(points, tangents, screw)?;
1023            Ok((r_axes, s_axes, None))
1024        }
1025        RingFrame::Planar(normal) => {
1026            let mut r_axes = Vec::with_capacity(points.len());
1027            let mut s_axes = Vec::with_capacity(points.len());
1028            for tangent in tangents {
1029                // Project out any tangent component and re-normalize, so the frame is
1030                // orthonormal even where the fitted plane and the sampled tangent
1031                // disagree in the last bits.
1032                let r = normal
1033                    .sub(tangent.scale(normal.dot(*tangent)))
1034                    .normalized()
1035                    .map_err(|_| {
1036                        "sweepSolid: the closed path's plane normal is parallel to its tangent"
1037                            .to_string()
1038                    })?;
1039                s_axes.push(tangent.cross(r).normalized()?);
1040                r_axes.push(r);
1041            }
1042            if requested.twist != 0.0 {
1043                roll_over_the_lap(points, &mut r_axes, &mut s_axes, requested.twist)?;
1044            }
1045            Ok((
1046                r_axes,
1047                s_axes,
1048                Some(SweepClosure {
1049                    holonomy: 0.0,
1050                    requested_twist: requested.twist,
1051                    applied_twist: requested.twist,
1052                    curve_shift: requested.shift,
1053                }),
1054            ))
1055        }
1056        RingFrame::Spatial => {
1057            let (r_axes, s_axes, closure) =
1058                ring_closure_frames(points, tangents, screw, requested)?;
1059            Ok((r_axes, s_axes, Some(closure)))
1060        }
1061    }
1062}
1063
1064/// The frames of a closed SPATIAL ring: the transported frame, counter-twisted
1065/// by its own holonomy so the section meets itself (see [`SweepClosure`]).
1066///
1067/// MEASURED, not assumed. The run is transported ONE STEP FURTHER than an open
1068/// run goes — from the last station back to station 0, by the very rule every
1069/// other step takes (the RMF, or a helix segment's screw frame under `Rigid`) —
1070/// and the returned axis is compared with station 0's: the signed angle between
1071/// them, right-handed about the start tangent, is the holonomy `h`.
1072///
1073/// Then station `k`, at chord-length `s_k` round a lap of chord length `L`
1074/// (the wrap chord included), is rolled about its tangent by `φ_k = −h·s_k/L`.
1075/// Rolling about the tangent commutes with the transport step (each step is a
1076/// rotation carrying one tangent onto the next, and it carries a roll about the
1077/// first onto the same roll about the second), so the frame transported from
1078/// the last twisted station across the wrap arrives rolled by
1079/// `h − h·s_{n−1}/L`, and station 0 sits `−h·(L − s_{n−1})/L` further on — the
1080/// same rate as every other span. The ring closes with no seam to skin.
1081///
1082/// No new station law: the closed lane already floors its budget at 16 stations
1083/// per quarter turn of the path's TURNING, and a closed curve turns at least
1084/// `2π` (Fenchel), so a lap carries at least 64 stations. The largest roll this
1085/// lays down is `π`, under 3° a station — finer than the 16 per quarter turn a
1086/// requested twist gets.
1087///
1088/// A REQUESTED twist `T` that closes ([`closing_twist`]) rides on the same rate:
1089/// each station rolls by `(T − h)·s_k/L`, so the section comes back from the lap
1090/// rolled by `T` alone, onto itself. [`sweep_profile_along_chain`] sizes the
1091/// station budget for `|T| + π`.
1092fn ring_closure_frames(
1093    points: &[Vec3],
1094    tangents: &[Vec3],
1095    screw: Option<(&[usize], &SweepPath)>,
1096    requested: ClosingTwist,
1097) -> Result<(Vec<Vec3>, Vec<Vec3>, SweepClosure), String> {
1098    let stations = points.len();
1099    if stations < 2 {
1100        return Err("sweepSolid: need at least 2 stations".into());
1101    }
1102    // The lap with its wrap: station 0 again, arriving as the LAST segment's end.
1103    let mut lap_points = points.to_vec();
1104    lap_points.push(points[0]);
1105    let mut lap_tangents = tangents.to_vec();
1106    lap_tangents.push(tangents[0]);
1107    let lap_segments: Option<Vec<usize>> = screw.map(|(segments, _)| {
1108        let mut extended = segments.to_vec();
1109        extended.push(segments[stations - 1]);
1110        extended
1111    });
1112    let lap_screw = match (&lap_segments, screw) {
1113        (Some(segments), Some((_, path))) => Some((segments.as_slice(), path)),
1114        _ => None,
1115    };
1116    let (mut r_axes, mut s_axes) = transported_frames(&lap_points, &lap_tangents, lap_screw)?;
1117    let returned = r_axes.pop().expect("a lap has a wrap station");
1118    s_axes.pop();
1119    let holonomy = returned.dot(s_axes[0]).atan2(returned.dot(r_axes[0]));
1120
1121    let applied = if requested.twist == 0.0 {
1122        -holonomy
1123    } else {
1124        requested.twist - holonomy
1125    };
1126    roll_over_the_lap(points, &mut r_axes, &mut s_axes, applied)?;
1127    Ok((
1128        r_axes,
1129        s_axes,
1130        SweepClosure {
1131            holonomy,
1132            requested_twist: requested.twist,
1133            applied_twist: applied,
1134            curve_shift: requested.shift,
1135        },
1136    ))
1137}
1138
1139/// Roll a ring's station frames about their tangents by `roll` over the lap,
1140/// linearly in chord length with the WRAP chord included: station `k` at `s_k`
1141/// of a lap of length `L` rolls by `roll·s_k/L`.
1142fn roll_over_the_lap(
1143    points: &[Vec3],
1144    r_axes: &mut [Vec3],
1145    s_axes: &mut [Vec3],
1146    roll: f64,
1147) -> Result<(), String> {
1148    let stations = points.len();
1149    let mut travelled = Vec::with_capacity(stations);
1150    let mut length = 0.0;
1151    for index in 0..stations {
1152        travelled.push(length);
1153        length += points[(index + 1) % stations].sub(points[index]).length();
1154    }
1155    if length <= 0.0 {
1156        return Err("sweepSolid: the closed path has zero length".into());
1157    }
1158    for station in 0..stations {
1159        let phi = roll * travelled[station] / length;
1160        let (sin_phi, cos_phi) = phi.sin_cos();
1161        let (r, s) = (r_axes[station], s_axes[station]);
1162        r_axes[station] = r.scale(cos_phi).add(s.scale(sin_phi));
1163        s_axes[station] = s.scale(cos_phi).sub(r.scale(sin_phi));
1164    }
1165    Ok(())
1166}
1167
1168/// The TRANSPORTED frames of a station run: seeded on an arbitrary perpendicular
1169/// to the first tangent and carried station to station — by the RMF, or by a
1170/// helix segment's screw frame under `Rigid` (see [`station_frames`]).
1171fn transported_frames(
1172    points: &[Vec3],
1173    tangents: &[Vec3],
1174    screw: Option<(&[usize], &SweepPath)>,
1175) -> Result<(Vec<Vec3>, Vec<Vec3>), String> {
1176    let stations = points.len();
1177    let mut r_axes = Vec::with_capacity(stations);
1178    let mut s_axes = Vec::with_capacity(stations);
1179    {
1180        let r0 = tangents[0].perpendicular()?; // any unit vector ⟂ T0
1181        s_axes.push(tangents[0].cross(r0).normalized()?);
1182        r_axes.push(r0);
1183        // The helix segment the running roll was fixed for, and that roll.
1184        let mut roll: Option<(usize, f64)> = None;
1185        for index in 0..stations - 1 {
1186            let t_next = tangents[index + 1];
1187            let helix = screw.and_then(|(segments, path)| {
1188                let segment = segments[index + 1];
1189                path.screw_axes[segment].map(|axis| (segment, axis, path))
1190            });
1191            let r_next = match helix {
1192                None => rmf_step(
1193                    r_axes[index],
1194                    tangents[index],
1195                    t_next,
1196                    points[index + 1].sub(points[index]),
1197                )
1198                .ok_or_else(|| format!("sweepSolid: frame degenerated at station {index}"))?,
1199                Some((segment, axis, path)) => {
1200                    // The roll is fixed ONCE per helix segment, at the station it
1201                    // departs from — station 0, or the joint the previous segment
1202                    // arrived at — so the incoming frame continues without a snap.
1203                    let phi = match roll {
1204                        Some((fixed, phi)) if fixed == segment => phi,
1205                        _ => {
1206                            let tangent = tangents[index];
1207                            let reference =
1208                                screw_reference(axis, tangent, path.name(segment), index)?;
1209                            let r = r_axes[index];
1210                            let phi = tangent.dot(reference.cross(r)).atan2(reference.dot(r));
1211                            roll = Some((segment, phi));
1212                            phi
1213                        }
1214                    };
1215                    let reference = screw_reference(axis, t_next, path.name(segment), index + 1)?;
1216                    let (sin_phi, cos_phi) = phi.sin_cos();
1217                    reference
1218                        .scale(cos_phi)
1219                        .add(t_next.cross(reference).scale(sin_phi))
1220                }
1221            };
1222            s_axes.push(t_next.cross(r_next).normalized()?);
1223            r_axes.push(r_next);
1224        }
1225    }
1226    Ok((r_axes, s_axes))
1227}
1228
1229/// A HELIX station's reference axis: the screw axis with its tangent component
1230/// removed, `normalize(W − (W·T)T)`.
1231///
1232/// Why this and not the transported frame: a helix is carried into itself by
1233/// ONE rigid motion — the screw about its axis — and that motion fixes both the
1234/// axis and the tangent's angle to it, so the axis' projection onto the normal
1235/// plane is the same body-fixed direction at every station. A frame built on it
1236/// IS the screw, so `R_k = F_k·F_0ᵀ` carries the profile about the coil's axis
1237/// and up its pitch. A rotation-minimizing frame has no spin about the tangent,
1238/// while the screw's frame spins at the helix's torsion `τ`, so the RMF rolls
1239/// away from it by `∫τ ds` — the rotation the user saw. On a TAPERED helix the
1240/// tangent's angle to the axis drifts with the radius, and this frame keeps the
1241/// section square to the axis all the same, which is what a tapered coil is.
1242///
1243/// Undefined only where the path runs ALONG its own axis, which a helix of
1244/// positive radius never does; refused by name there rather than guessed.
1245fn screw_reference(axis: Vec3, tangent: Vec3, segment: &str, station: usize) -> Result<Vec3, String> {
1246    axis.sub(tangent.scale(axis.dot(tangent)))
1247        .normalized()
1248        .map_err(|_| {
1249            format!(
1250                "sweepSolid: path segment '{segment}' runs along its own helix axis at station \
1251                 {station}, so the screw motion that carries the section has no direction \
1252                 about the path there"
1253            )
1254        })
1255}
1256
1257/// The stations a path is sampled at, and what each one carries: the point, the
1258/// unit tangent, and the CURVATURE VECTOR `κ·N̂` — the path's own bend there,
1259/// measured from the curve rather than differenced off the station polyline,
1260/// because it is the number the tight-bend refusal quotes.
1261pub(crate) struct StationSamples {
1262    pub(crate) points: Vec<Vec3>,
1263    pub(crate) tangents: Vec<Vec3>,
1264    /// `κ·N̂` per station: its length is the curvature and its direction points
1265    /// at the centre of the osculating circle. Zero on a straight piece.
1266    pub(crate) curvatures: Vec<Vec3>,
1267    /// The path SEGMENT each station belongs to. A joint station is emitted by
1268    /// the segment ARRIVING at it, so it carries that segment's index.
1269    pub(crate) segments: Vec<usize>,
1270}
1271
1272/// The CURVATURE VECTOR of a curve from its first two derivatives:
1273///
1274/// ```text
1275///     κ·N̂ = (r̈ − (r̈·T̂)T̂) / |ṙ|²,      T̂ = ṙ/|ṙ|
1276/// ```
1277///
1278/// the component of the second derivative perpendicular to the tangent over the
1279/// squared speed — the reparameterization-invariant form, so a non-uniformly
1280/// parameterized path reads its own geometric bend rather than its
1281/// parameterization's. A degree-1 piece has a zero second derivative
1282/// (`derivatives` pads past the degree), so a straight segment reads exactly
1283/// zero and needs no special case.
1284pub(crate) fn curvature_vector(first: Vec3, second: Vec3) -> Vec3 {
1285    let speed_squared = first.dot(first);
1286    if speed_squared <= 0.0 {
1287        return Vec3::default();
1288    }
1289    let tangent = first.scale(1.0 / speed_squared.sqrt());
1290    second
1291        .sub(tangent.scale(second.dot(tangent)))
1292        .scale(1.0 / speed_squared)
1293}
1294
1295/// Sample a CHAIN as ONE trajectory: the per-segment arc-length estimate, the
1296/// station budget it implies, and the stations themselves.
1297///
1298/// Extracted from the chain sweep so the tight-bend INSTRUMENT
1299/// ([`sweep_bend_profile`]) measures the very stations the sweep lofts through,
1300/// rather than a second sampling that agrees with it today.
1301fn chain_stations(path: &SweepPath, stations: usize) -> Result<StationSamples, String> {
1302    let tolerance = 1e-6;
1303    let chain = path.curves.as_slice();
1304    // --- Per-segment arc-length estimate, and the domain each one is sampled
1305    //     over. A degenerate domain is refused here rather than producing a
1306    //     silently collapsed station block.
1307    //
1308    //     Each segment's TURNING — the total angle its tangent sweeps — is measured
1309    //     too, because that, not its length, is what decides how many stations it
1310    //     needs (see the floor below).
1311    let mut domains = Vec::with_capacity(chain.len());
1312    let mut lengths = Vec::with_capacity(chain.len());
1313    let mut turning = Vec::with_capacity(chain.len());
1314    for (index, curve) in chain.iter().enumerate() {
1315        let [t0, t1] = curve.domain()?;
1316        if (t1 - t0).abs() <= tolerance {
1317            return Err(format!(
1318                "sweepSolid: path segment {index} has a degenerate domain"
1319            ));
1320        }
1321        let mut length = 0.0;
1322        let mut previous = curve.evaluate(t0)?;
1323        for sample in 1..=16 {
1324            let point = curve.evaluate(t0 + (t1 - t0) * sample as f64 / 16.0)?;
1325            length += point.sub(previous).length();
1326            previous = point;
1327        }
1328        domains.push([t0, t1]);
1329        lengths.push(length);
1330        turning.push(segment_turning(curve, [t0, t1]).map_err(|error| {
1331            format!("sweepSolid: path tangent is degenerate on segment {index}: {error}")
1332        })?);
1333    }
1334    let total_length: f64 = lengths.iter().sum();
1335    if total_length <= tolerance {
1336        return Err("sweepSolid: path chain has zero length".into());
1337    }
1338
1339    // --- Station budget per segment: proportional to arc length, floored at 2
1340    //     (a segment needs its two ends), and the remainder handed to the
1341    //     longest segment so the totals land exactly on `stations`.
1342    //
1343    //     The floor is raised to the segment's own TURNING, at 16
1344    //     stations per quarter turn — the density the twisted builder already
1345    //     uses, and the one that puts the cubic v-interpolant's own error around
1346    //     1e-7 relative (it falls as the fourth power of the angular spacing,
1347    //     `Δθ⁴/384`, which is ~2.4e-7 at 16 per quarter turn).
1348    //
1349    //     WHY THE CLOSED LANE NEEDS IT AND THE OPEN LANE DOES NOT, stated as the
1350    //     measurement that forced it: a rounded rectangle's four quarter arcs are
1351    //     SHORT next to its four straight runs, so a length-proportional share of
1352    //     32 hands each arc the 2-station floor — its two ends and nothing in
1353    //     between — and the interpolant cuts the corner. Measured on a 30 x 18
1354    //     path with r = 4: the ring came out 6.2 % under `A · perimeter`, a
1355    //     watertight, validating, WRONG solid. The floor takes the same case to
1356    //     1e-5.
1357    //
1358    //     An OPEN path is under-sampled the same way and now takes the same floor.
1359    //     It did not until 2026-09-25, and a coil showed what that cost: an 8-turn
1360    //     helix got 32 stations, four a turn, and the loft between sections a
1361    //     quarter turn apart folded its own wall — a `Rigid` sweep along the coil,
1362    //     or a coil with a straight lead, was refused outright. (A curvature-JUMP
1363    //     law at G1 joints is still this plan's remaining item 2; the floor is
1364    //     about how far a segment turns, not how its curvature changes.)
1365    //
1366    //     Once the floors lift the turning segments past their share, the SHARES
1367    //     are taken of the lifted total, not of the `stations` asked for: the
1368    //     segments that do not turn are sampled at the density the turning ones
1369    //     were raised to. Measured on the coil above with a tangent straight lead:
1370    //     the floor gave the coil ~500 stations and the lead its share of 32 — its
1371    //     two ends — and the interpolant, stiff from the coil, bowed through the
1372    //     lead: 2.7e-5 off Pappus on a lead of 10, 5.0e-4 on a lead of 100, and
1373    //     1.5e-7 / 2.1e-7 at the density the coil was given.
1374    let floors: Vec<usize> = turning
1375        .iter()
1376        .map(|turn| ((turn / std::f64::consts::FRAC_PI_2) * 16.0).ceil() as usize)
1377        .collect();
1378    let floored_length: f64 = lengths
1379        .iter()
1380        .zip(&floors)
1381        .filter(|(_, floor)| **floor > 2)
1382        .map(|(length, _)| length)
1383        .sum();
1384    let floor_total: usize = floors.iter().sum();
1385    //     A CLOSED path is not lifted — its shares stay those of the `stations`
1386    //     asked for (its floors can still move where a segment's turning is
1387    //     re-read densely, see `segment_turning`): its loft
1388    //     runs every skin through as many laps as its section curves take to
1389    //     cycle, so a lifted total multiplies against the cap there, and the
1390    //     closed lane's measured residuals are pinned at its own budget.
1391    let lifted = if !path.closed && floored_length > tolerance {
1392        let straight = floor_total as f64 * (total_length - floored_length) / floored_length;
1393        ((floor_total as f64 + straight).ceil() as usize).min(1024)
1394    } else {
1395        0
1396    };
1397    let stations = stations.max(lifted);
1398    let mut budget: Vec<usize> = lengths
1399        .iter()
1400        .zip(&floors)
1401        .map(|(length, floor)| {
1402            let share = (stations as f64 * length / total_length).round() as usize;
1403            share.max(*floor).max(2)
1404        })
1405        .collect();
1406    // Each interior joint shares one station between its two segments, so the
1407    // emitted total is sum(budget) - (segments - 1). A CLOSED path drops its very
1408    // last station too — it is the first one — so it emits one fewer again.
1409    let shared = chain.len() - 1 + usize::from(path.closed);
1410    let emitted: usize = budget.iter().sum::<usize>() - shared;
1411    // The 2-per-segment floor means a chain of n segments emits at least n+1
1412    // stations, so a long enough chain cannot be held under the loft's station
1413    // cap by any distribution. The loft's interpolation solve is O(stations³)
1414    // per control column, so silently accepting such a path would look like a
1415    // hang. Refuse with the count instead — the same honesty the twist cap uses.
1416    const MAX_STATIONS: usize = 1024;
1417    if emitted > MAX_STATIONS {
1418        return Err(format!(
1419            "sweepSolid: a {}-segment path needs at least {emitted} stations, past the \
1420             {MAX_STATIONS}-station cap the loft solve can carry; sweep it in fewer, \
1421             longer pieces",
1422            chain.len()
1423        ));
1424    }
1425    if emitted < stations {
1426        let longest = lengths
1427            .iter()
1428            .enumerate()
1429            .max_by(|a, b| a.1.total_cmp(b.1))
1430            .map(|(index, _)| index)
1431            .unwrap_or(0);
1432        budget[longest] += stations - emitted;
1433    }
1434
1435    // --- Sample the chain as ONE trajectory.
1436    let mut points: Vec<Vec3> = Vec::new();
1437    let mut tangents: Vec<Vec3> = Vec::new();
1438    let mut curvatures: Vec<Vec3> = Vec::new();
1439    let mut segments: Vec<usize> = Vec::new();
1440    for (index, curve) in chain.iter().enumerate() {
1441        let [t0, t1] = domains[index];
1442        let count = budget[index];
1443        // Skip the first sample of every segment after the first: that station
1444        // is the joint, already emitted as the previous segment's last. On a
1445        // CLOSED path the LAST segment's last sample is skipped too — it is
1446        // station 0, and a ring that emitted it twice would interpolate a
1447        // zero-length span across its own seam.
1448        let first = usize::from(index > 0);
1449        if first == 1 {
1450            // THE JOINT reads the TIGHTER of the two pieces that meet there. A
1451            // G1 joint is tangent-continuous, not curvature-continuous — a line
1452            // running into a cubic hands over at κ = 0 and the cubic starts at
1453            // its own spike — and the station is emitted once, by the piece that
1454            // ARRIVES. Left as that, a path whose tightest bend sits at the
1455            // START of a curved piece would quote the straight piece's κ = 0 and
1456            // the bend would go unnamed; on the wire harness's own span the two
1457            // ends are mirror images and only ONE of them was being read.
1458            let derivatives = curve.derivatives_small(t0, 2)?;
1459            let arriving = curvature_vector(derivatives[1], derivatives[2]);
1460            if let Some(departing) = curvatures.last_mut() {
1461                if arriving.length() > departing.length() {
1462                    *departing = arriving;
1463                }
1464            }
1465        }
1466        let last = count - usize::from(path.closed && index + 1 == chain.len());
1467        for sample in first..last {
1468            let t = t0 + (t1 - t0) * sample as f64 / (count - 1) as f64;
1469            let derivatives = curve.derivatives_small(t, 2)?;
1470            let tangent = derivatives[1].normalized().map_err(|_| {
1471                format!("sweepSolid: path tangent is degenerate on segment {index}")
1472            })?;
1473            points.push(derivatives[0]);
1474            tangents.push(tangent);
1475            curvatures.push(curvature_vector(derivatives[1], derivatives[2]));
1476            segments.push(index);
1477        }
1478    }
1479    // A CLOSED path's seam is a joint like any other, and station 0 is the only
1480    // one whose two sides are sampled at opposite ends of the loop: the last
1481    // segment's end is station 0, and it was dropped rather than emitted.
1482    if path.closed {
1483        if let (Some(curve), Some(domain)) = (chain.last(), domains.last()) {
1484            let derivatives = curve.derivatives_small(domain[1], 2)?;
1485            let arriving = curvature_vector(derivatives[1], derivatives[2]);
1486            if let Some(departing) = curvatures.first_mut() {
1487                if arriving.length() > departing.length() {
1488                    *departing = arriving;
1489                }
1490            }
1491        }
1492    }
1493    Ok(StationSamples {
1494        points,
1495        tangents,
1496        curvatures,
1497        segments,
1498    })
1499}
1500
1501/// The ring frame's fixed reference axis: the path plane's own normal, oriented
1502/// once and reused at every station.
1503///
1504/// A planar path's rotation-minimizing frame IS the planar frame — the transport
1505/// has no component about the tangent, so the axis perpendicular to the plane is
1506/// carried unchanged the whole way round and comes back to itself EXACTLY. That
1507/// is what makes a closed planar sweep need no twist policy: the holonomy is
1508/// zero by construction, not by cancellation. Taking the plane normal as the
1509/// reference axis says so in the construction rather than hoping a
1510/// station-by-station transport returns to where it started.
1511///
1512/// The guard is the one case a fixed axis cannot serve: a tangent parallel to the
1513/// normal, which a planar path cannot have and a mis-fitted plane could.
1514fn ring_normal(normal: Vec3, tangents: &[Vec3]) -> Result<Vec3, String> {
1515    for (index, tangent) in tangents.iter().enumerate() {
1516        if tangent.cross(normal).length() < 0.5 {
1517            return Err(format!(
1518                "sweepSolid: the closed path's own plane normal is within 30° of its tangent \
1519                 at station {index}, so the path does not lie in that plane after all"
1520            ));
1521        }
1522    }
1523    normal.normalized().map_err(|_| {
1524        "sweepSolid: the closed path's plane normal is degenerate".to_string()
1525    })
1526}
1527
1528/// The total TURNING of one path segment: the angle its unit tangent sweeps over
1529/// its own domain, summed over 16 uniform steps — or, where a step turns more
1530/// than [`TURNING_ALIAS_STEP`], over 16 steps per knot span.
1531///
1532/// Not the same question as its length, and the one a station budget has to
1533/// answer: a straight segment of any length is carried exactly by two stations,
1534/// while a quarter arc of any radius needs enough stations that the cubic
1535/// interpolant through them does not cut the corner. 16 samples is the density
1536/// every other estimate in this builder uses — and it ALIASES a segment that
1537/// turns more than a few times: on an 8-turn coil each step is half a turn, and
1538/// the angle between two tangents cannot read more than π however far the curve
1539/// turned between them. A step that large is the sign, so only then is the
1540/// segment re-read densely; everything that reads calmly at 16 keeps exactly the
1541/// reading, and so exactly the floor, it always had (a half circle reads a hair
1542/// over π there, and a floor of 33 — the budget the closed lane is pinned at).
1543fn segment_turning(curve: &NurbsCurve, domain: [f64; 2]) -> Result<f64, String> {
1544    let (total, widest) = turning_at(curve, domain, 16)?;
1545    if widest <= TURNING_ALIAS_STEP {
1546        return Ok(total);
1547    }
1548    let spans = curve.knots.windows(2).filter(|pair| pair[1] > pair[0]).count();
1549    Ok(turning_at(curve, domain, (spans * 16).clamp(64, 8192))?.0)
1550}
1551
1552/// The widest tangent step [`segment_turning`] trusts at 16 samples: π/8, twice
1553/// the step of a quarter arc read at that density.
1554const TURNING_ALIAS_STEP: f64 = std::f64::consts::PI / 8.0;
1555
1556/// `(total turning, widest single step)` over `steps` uniform steps.
1557fn turning_at(curve: &NurbsCurve, domain: [f64; 2], steps: usize) -> Result<(f64, f64), String> {
1558    let [t0, t1] = domain;
1559    let mut total = 0.0;
1560    let mut widest: f64 = 0.0;
1561    // A stationary sample — a cusp in the parameterization, which an involute
1562    // flank's Hermite fit carries by construction — still has a direction of
1563    // travel, and the turning is a question about directions.
1564    let mut previous = curve.unit_tangent(t0, t0, t1)?;
1565    for sample in 1..=steps {
1566        let tangent = curve.unit_tangent(t0 + (t1 - t0) * sample as f64 / steps as f64, t0, t1)?;
1567        let step = previous.dot(tangent).clamp(-1.0, 1.0).acos();
1568        total += step;
1569        widest = widest.max(step);
1570        previous = tangent;
1571    }
1572    Ok((total, widest))
1573}
1574
1575/// ONE double-reflection step (Wang et al. 2008) of a rotation-minimizing frame:
1576/// carry the reference axis `r` from a station with tangent `t` to the next with
1577/// tangent `t_next` across the step vector `step`, then re-orthogonalize against
1578/// the new tangent to shed floating drift. `None` when the result degenerates.
1579///
1580/// Extracted so the frame the sections are placed on and the frame the holonomy
1581/// is measured with are the SAME transport rather than two copies of it.
1582fn rmf_step(r: Vec3, t: Vec3, t_next: Vec3, step: Vec3) -> Option<Vec3> {
1583    let c1 = step.dot(step);
1584    let candidate = if c1 <= 1e-18 {
1585        // Coincident stations: carry the reference axis forward unchanged.
1586        r
1587    } else {
1588        // First reflection across the plane bisecting the step vector.
1589        let reflected_r = r.sub(step.scale(2.0 / c1 * step.dot(r)));
1590        let reflected_t = t.sub(step.scale(2.0 / c1 * step.dot(t)));
1591        // Second reflection across the plane bisecting the tangents.
1592        let v2 = t_next.sub(reflected_t);
1593        let c2 = v2.dot(v2);
1594        if c2 <= 1e-18 {
1595            reflected_r
1596        } else {
1597            reflected_r.sub(v2.scale(2.0 / c2 * v2.dot(reflected_r)))
1598        }
1599    };
1600    candidate
1601        .sub(t_next.scale(candidate.dot(t_next)))
1602        .normalized()
1603        .ok()
1604}
1605
1606// ===========================================================================
1607// The TIGHT BEND: a section carried through a bend tighter than its own reach
1608// ===========================================================================
1609
1610/// How many points each profile curve contributes to the section's REACH.
1611///
1612/// The reach has to be a WITNESS — a point of the section that really does turn
1613/// back — so it is read off the curve itself and never off the control hull:
1614/// `make_arc` places a half-circle's middle controls at `r√2`, so a hull bound
1615/// would refuse a round section 41 % early. Sampling under-reads instead, which
1616/// is the safe direction: at 64 samples per curve a circle's support is short by
1617/// `r(1 − cos(π/128))` = 3.0e-4·r, so the guard can only ever refuse LATE, and a
1618/// bend that slips past it still meets the fold acceptance at the loft's exit.
1619const BEND_SECTION_SAMPLES: usize = 64;
1620
1621/// One station's tight-bend reading: what the path does there and how far the
1622/// section reaches into it.
1623///
1624/// The quantity that decides the sweep is [`BendStation::ratio`] — `reach·κ`,
1625/// the same `ρκ` the fillet's canal surface folds at. A section point offset `d`
1626/// towards the centre of the osculating circle advances `(1 − d·κ)` per unit of
1627/// path, so at `d·κ = 1` it stops and past it it goes BACKWARDS: the swept
1628/// envelope has turned back through itself and the solid is the envelope trimmed
1629/// at that self-intersection, not the loft through the crossing sections.
1630#[derive(Debug, Clone, Copy)]
1631pub struct BendStation {
1632    /// The station's index in the sweep's own station list.
1633    pub station: usize,
1634    /// Where it is on the path.
1635    pub point: Vec3,
1636    /// The path's curvature `κ` there (zero on a straight piece).
1637    pub curvature: f64,
1638    /// The path's radius of curvature `ρ = 1/κ`, `f64::INFINITY` where straight.
1639    pub radius: f64,
1640    /// How far the placed section reaches towards the centre of the bend —
1641    /// its support in the curvature-normal direction. Negative when the whole
1642    /// section sits on the OUTSIDE of the bend, zero where the path is straight.
1643    pub reach: f64,
1644    /// `reach·κ`. `< 1` builds; `≥ 1` folds.
1645    pub ratio: f64,
1646}
1647
1648/// The profile's sample points in the placement frame's own 2D coordinates
1649/// `(a, b) = ((x − origin)·pu, (x − origin)·pv)`.
1650///
1651/// Sampled ONCE: `Transplant` places the same 2D shape at every station
1652/// (`world = station + rᵢ·a + sᵢ·b`), so a station's reach is a pair of dot
1653/// products against the precomputed cloud rather than a re-evaluation of the
1654/// placed curves.
1655fn section_offsets(
1656    profile: &[NurbsCurve],
1657    placement: &ProfileAnchor,
1658) -> Result<Vec<(f64, f64)>, String> {
1659    let ProfileAnchor { origin, pu, pv, .. } = *placement;
1660    let mut offsets = Vec::with_capacity(profile.len() * BEND_SECTION_SAMPLES);
1661    for curve in profile {
1662        let [start, end] = curve.domain()?;
1663        for sample in 0..BEND_SECTION_SAMPLES {
1664            let t = start + (end - start) * sample as f64 / BEND_SECTION_SAMPLES as f64;
1665            let local = curve.evaluate(t)?.sub(origin);
1666            offsets.push((local.dot(pu), local.dot(pv)));
1667        }
1668    }
1669    Ok(offsets)
1670}
1671
1672/// The `ρκ` profile: one [`BendStation`] per station, in station order.
1673///
1674/// `axes` are the frames the sections are PLACED on (twist included), because
1675/// the reach is the section's support in the curvature-normal direction and a
1676/// non-round section rolled about the path presents a different support there.
1677fn bend_stations(
1678    points: &[Vec3],
1679    curvatures: &[Vec3],
1680    axes: &[(Vec3, Vec3)],
1681    offsets: &[(f64, f64)],
1682) -> Vec<BendStation> {
1683    (0..points.len())
1684        .map(|station| {
1685            let curvature = curvatures[station].length();
1686            if curvature <= 0.0 {
1687                return BendStation {
1688                    station,
1689                    point: points[station],
1690                    curvature: 0.0,
1691                    radius: f64::INFINITY,
1692                    reach: 0.0,
1693                    ratio: 0.0,
1694                };
1695            }
1696            // The curvature vector is perpendicular to the tangent by
1697            // construction, so its components on (r, s) are the whole of it.
1698            let normal = curvatures[station].scale(1.0 / curvature);
1699            let (r_axis, s_axis) = axes[station];
1700            let (nu, nv) = (normal.dot(r_axis), normal.dot(s_axis));
1701            let reach = offsets
1702                .iter()
1703                .map(|(a, b)| a * nu + b * nv)
1704                .fold(f64::NEG_INFINITY, f64::max);
1705            BendStation {
1706                station,
1707                point: points[station],
1708                curvature,
1709                radius: 1.0 / curvature,
1710                reach,
1711                ratio: reach * curvature,
1712            }
1713        })
1714        .collect()
1715}
1716
1717/// The `ρκ` profile a TRANSPLANT sweep of `profile` along `path` would carry at
1718/// `stations` stations — the instrument behind the refusal below, exposed so a
1719/// caller (or a probe) can measure a path BEFORE asking for the solid.
1720///
1721/// Reads the very stations the sweep lofts through: the same budget, the same
1722/// sampling, the same rotation-minimizing frame. Untwisted, which is the lane
1723/// every caller of [`sweep_profile_along_chain_with_stations`] is in.
1724pub fn sweep_bend_profile(
1725    profile: &[NurbsCurve],
1726    path: &SweepPath,
1727    stations: usize,
1728) -> Result<Vec<BendStation>, String> {
1729    let placement = resolve_placement(profile, None, SectionPlacement::Transplant)?;
1730    let StationSamples {
1731        points,
1732        tangents,
1733        curvatures,
1734        ..
1735    } = chain_stations(path, stations.clamp(2, 1024))?;
1736    let closure = ring_frame(path, &tangents)?;
1737    let untwisted = ClosingTwist {
1738        twist: 0.0,
1739        shift: 0,
1740    };
1741    let (r_axes, s_axes, _) = station_frames(&points, &tangents, closure, None, untwisted)?;
1742    let axes: Vec<(Vec3, Vec3)> = r_axes.into_iter().zip(s_axes).collect();
1743    let offsets = section_offsets(profile, &placement)?;
1744    Ok(bend_stations(&points, &curvatures, &axes, &offsets))
1745}
1746
1747/// The prefix every tight-bend refusal starts with, so a caller can CLASSIFY one
1748/// (the wire harness reports it as a bundle status) without matching free text.
1749pub const SWEEP_TIGHT_BEND_REFUSAL: &str = "sweepSolid: the path bends TIGHTER than the section";
1750
1751/// The refusal for a path that bends tighter than the section reaches, or `None`
1752/// when every station builds.
1753///
1754/// Quotes the WORST station rather than the first one scanned — the tightest
1755/// bend is the one the user has to ease, and a path can cross the bar in several
1756/// places (both ends of one Hermite span, typically).
1757fn tight_bend_refusal(profile: &[BendStation], stations: usize) -> Option<String> {
1758    let folded = profile.iter().filter(|station| station.ratio >= 1.0).count();
1759    if folded == 0 {
1760        return None;
1761    }
1762    let worst = profile
1763        .iter()
1764        .max_by(|a, b| a.ratio.total_cmp(&b.ratio))
1765        .expect("a folded station exists");
1766    Some(format!(
1767        "{SWEEP_TIGHT_BEND_REFUSAL} at station {} of {stations} ({:.4}, {:.4}, {:.4}): the path's \
1768         radius of curvature there is {:.6} and the section reaches {:.6} into the bend, so \
1769         ρ·κ = {:.3} ≥ 1 ({folded} of {stations} stations do). A section carried through a bend \
1770         tighter than its own reach sweeps an envelope that passes through ITSELF — the solid \
1771         there is that envelope trimmed at its own self-intersection, which this lane does not \
1772         build: it skins the section between stations, and past ρ·κ = 1 the sections CROSS and \
1773         the wall folds. Ease the bend to a radius above {:.6}, or narrow the section until it \
1774         reaches no further than {:.6} (a round section of diameter {:.6})",
1775        worst.station,
1776        worst.point.x,
1777        worst.point.y,
1778        worst.point.z,
1779        worst.radius,
1780        worst.reach,
1781        worst.ratio,
1782        worst.reach,
1783        worst.radius,
1784        2.0 * worst.radius,
1785    ))
1786}
1787
1788/// The part of the sweep that does not care HOW the stations were sampled:
1789/// validate the profile, distribute the twist over the sampled polyline,
1790/// propagate the rotation-minimizing frame, place a profile copy per station
1791/// and loft through them.
1792///
1793/// Split out so a MULTI-SEGMENT path can share it
1794/// ([`sweep_profile_along_chain_stations`]). The single-curve sampler above is
1795/// deliberately left calling this with its OWN sampling rather than being
1796/// re-expressed as a one-element chain: two code paths that agree today are not
1797/// the same thing as one code path, and the untwisted single-curve geometry is
1798/// pinned by volume tests that must not move.
1799///
1800/// The RMF is propagated over the station sequence as a whole, so a chain's
1801/// frame carries across a joint exactly as it carries across any other pair of
1802/// adjacent stations — there is no per-segment restart to reconcile.
1803///
1804/// `closure` turns the run into a RING: anything but [`RingFrame::Open`] says the
1805/// stations come back round to station 0 (which is therefore NOT repeated at the
1806/// end). `Planar(normal)` says the path lies in the plane of that normal, so the
1807/// reference axis is fixed rather than transported — a planar frame, whose
1808/// holonomy round the loop is zero by construction. `Spatial` transports the
1809/// frame and closes it with the counter-twist policy ([`SweepClosure`]). Either
1810/// way the sections skin through a CLOSED loft: one face per section curve, no
1811/// end caps, `V − E + F = 0`. `Open` is the open run: transported frame, two end
1812/// caps, and the refusal below for a trajectory that closes anyway.
1813///
1814/// `chain` is the station-to-segment map and the classified path, from the chain
1815/// sampler; the single-curve sampler has none. Only `Rigid` reads it, to frame a
1816/// HELIX segment's stations off its screw axis ([`screw_reference`]).
1817#[allow(clippy::too_many_arguments)]
1818fn sweep_sections_through_samples(
1819    profile: &[NurbsCurve],
1820    points: &[Vec3],
1821    tangents: &[Vec3],
1822    curvatures: &[Vec3],
1823    twist_angle: f64,
1824    placement: ProfileAnchor,
1825    placement_mode: SectionPlacement,
1826    closure: RingFrame,
1827    envelope_gap: Option<&str>,
1828    chain: Option<(&[usize], &SweepPath)>,
1829) -> Result<BrepSolid, String> {
1830    let tolerance = 1e-6;
1831    let stations = points.len();
1832    if stations < 2 {
1833        return Err("sweepSolid: need at least 2 stations".into());
1834    }
1835
1836    let ProfileAnchor { origin, pu, pv, .. } = placement;
1837
1838    // --- 2b. Twist distribution: cumulative ARC-LENGTH fractions over the
1839    //         sampled station polyline (chord sums), so the twist advances
1840    //         uniformly in space even on a non-uniformly parameterized path.
1841    //         Only computed when a twist is actually requested — the
1842    //         `twist_angle == 0.0` path must stay bit-identical to the
1843    //         pre-twist builder — and only on an OPEN run. These fractions stop at
1844    //         the last station, so a RING laid down by them would roll its wrap
1845    //         span at a different rate from every other; a ring rolls its twist
1846    //         in its frames instead, over the whole lap (§3).
1847    let twist_fractions: Option<Vec<f64>> = if twist_angle != 0.0 && !closure.is_ring() {
1848        let mut cumulative = vec![0.0; stations];
1849        let mut total = 0.0;
1850        for index in 1..stations {
1851            total += points[index].sub(points[index - 1]).length();
1852            cumulative[index] = total;
1853        }
1854        if total <= tolerance {
1855            return Err("sweepSolid: path has zero length; cannot distribute the twist".into());
1856        }
1857        for length in &mut cumulative {
1858            *length /= total;
1859        }
1860        Some(cumulative)
1861    } else {
1862        None
1863    };
1864
1865    // --- 3. The station frames.
1866    //
1867    //   OPEN: rotation-minimizing frames via the double-reflection method, seeded
1868    //   on an arbitrary perpendicular to T0 — arbitrary because `Rigid` cancels
1869    //   the seed out (`R_k = F_k·F_0ᵀ`) and `Transplant` has always taken whatever
1870    //   roll `perpendicular()` returned.
1871    //
1872    //   RING: the plane normal is the reference axis AT EVERY STATION. A planar
1873    //   path's RMF is that frame already — the transport has no component about
1874    //   the tangent — so this is not a different frame from the open lane's, it is
1875    //   the same frame written down in closed form instead of accumulated. Which
1876    //   matters for exactly one reason: written down, it returns to itself at the
1877    //   seam EXACTLY, so the ring's last section meets its first with no residual
1878    //   roll to skin away.
1879    //
1880    //   SPATIAL RING: the open lane's transported frame, carried one step further
1881    //   round the wrap to MEASURE the roll it returns with (the holonomy), and
1882    //   counter-twisted by that roll linearly in arc length so the lap's last span
1883    //   rolls at the same rate as every other and the ring meets itself
1884    //   (`ring_closure_frames`).
1885    //
1886    //   A REQUESTED twist on a RING is laid down in those frames too, over the
1887    //   whole lap with the wrap, on top of the counter-twist — but only a twist
1888    //   that closes: after the lap the section comes back rolled by it, so it must
1889    //   carry the section onto itself ([`closing_twist`], refused by name
1890    //   otherwise). A fractional turn on a symmetric section closes onto a
1891    //   DIFFERENT curve than it left, `shift` along, and the loft below is told.
1892    //
1893    //   RIGID on a HELIX segment: the helix's own screw frame (see
1894    //   `screw_reference`), which the RMF rolls away from by the helix's
1895    //   integrated torsion. `Transplant` keeps the RMF everywhere — it squares
1896    //   the section to the path and centres it there, and its callers (SWP, the
1897    //   tube, the harness) are not this contract.
1898    let screw = chain.filter(|_| placement_mode == SectionPlacement::Rigid);
1899    let requested = if closure.is_ring() && twist_angle != 0.0 {
1900        let (section, _) = place_section(profile, placement, placement_mode, points[0], tangents[0])?;
1901        ring_closing_twist(&section, points, tangents, twist_angle)?
1902    } else {
1903        ClosingTwist {
1904            twist: 0.0,
1905            shift: 0,
1906        }
1907    };
1908    let (r_axes, s_axes, _) = station_frames(points, tangents, closure, screw, requested)?;
1909
1910    // --- 3b. An OPEN trajectory that comes back to where it started cannot be
1911    //         capped: the first and last sections land on top of each other, and
1912    //         the loft's own complaint about that ("end sections coincide") says
1913    //         nothing about the path. Name it here instead, with the two things
1914    //         that do work. The check is on the PATH's first and last station
1915    //         POINTS, which is enough: a path that comes back to where it started
1916    //         puts the two caps in the same place whatever direction it arrives
1917    //         from — coincident when it also arrives pointing the same way,
1918    //         interpenetrating when it does not, and neither is a solid worth
1919    //         building.
1920    //
1921    //         TWO CLOSURE NOTIONS COEXIST, and deliberately. A path whose ends
1922    //         meet within the chainer's join band is CLASSIFIED closed and swept
1923    //         as a ring above — that is the shape a user drew as a loop. This band
1924    //         is a thousandth of the distance TRAVELLED, a hundred times looser,
1925    //         and it catches the run that is nearly a loop without being one: a
1926    //         359.7° arc whose caps would interpenetrate rather than coincide. The
1927    //         first is a shape with a construction; the second is an accident, and
1928    //         it is still refused. On a circular arc the band is the last ~0.36°
1929    //         before closure, so a 359° sweep still builds.
1930    if !closure.is_ring() {
1931        let travel: f64 = points
1932            .windows(2)
1933            .map(|pair| pair[1].sub(pair[0]).length())
1934            .sum();
1935        let gap = points[stations - 1].sub(points[0]).length();
1936        if travel > tolerance && gap <= 1e-3 * travel {
1937            return Err(format!(
1938                "sweepSolid: the path returns to where it started ({gap:.3e} apart after \
1939                 travelling {travel:.3e}), so the sweep's two end caps would land on top of each \
1940                 other. A closed path swept this way IS a revolution — build it with Revolve \
1941                 about the same axis, or close the path exactly so it sweeps as a capless ring, \
1942                 or sweep the run in two halves and union them"
1943            ));
1944        }
1945    }
1946
1947    // --- 4. Place a transformed copy of the profile at each station, by the mode
1948    //        the caller asked for (see `SectionPlacement`).
1949    //
1950    //   TRANSPLANT maps the profile's own plane frame onto the station frame,
1951    //       world = P_k + r_k·(l·pu) + s_k·(l·pv),   l = x − origin
1952    //   which lands the anchor origin ON the path and the profile square to it.
1953    //   The profile's distance from the path and its angle to the path are both
1954    //   discarded — there is no term along the tangent to carry them.
1955    //
1956    //   RIGID applies the path's own motion to the profile where it stands,
1957    //       world = P_k + a·T_k + b·r_k + c·s_k,     (a,b,c) = (l·T_0, l·r_0, l·s_0)
1958    //                                                l = x − P_0
1959    //   which is `P_k + R_k·(x − P_0)` for `R_k = F_k·F_0ᵀ`, written in the frame
1960    //   basis rather than as a matrix. The tangent term `a·T_k` is exactly what
1961    //   TRANSPLANT drops, and carrying it is what preserves both the profile's
1962    //   offset from the path and its angle to it. At station 0 the frames
1963    //   coincide and the section IS the profile, so the START cap lies in the
1964    //   plane it was drawn in.
1965    let path_start = points[0];
1966    let (t0, r0, s0) = (tangents[0], r_axes[0], s_axes[0]);
1967    let rigid = placement_mode == SectionPlacement::Rigid;
1968    let mut sections: Vec<Vec<NurbsCurve>> = Vec::with_capacity(stations);
1969    // Only the RIGID guard below reads these; TRANSPLANT leaves them empty.
1970    let mut station_points: Vec<Vec<Vec3>> = Vec::with_capacity(if rigid { stations } else { 0 });
1971    let mut section_normals: Vec<Vec3> = Vec::with_capacity(if rigid { stations } else { 0 });
1972    let mut anchor_track: Vec<Vec3> = Vec::with_capacity(if rigid { stations } else { 0 });
1973    // The axes the section is actually PLACED on, twist included — what the
1974    // tight-bend guard below measures the section's reach in.
1975    let mut placed_axes: Vec<(Vec3, Vec3)> = Vec::with_capacity(stations);
1976    for station in 0..stations {
1977        let station_origin = points[station];
1978        // Rotate the RMF axes about the tangent by the station's twist angle
1979        // (Rodrigues on vectors ⟂ the tangent: r' = r·cosφ + s·sinφ,
1980        // s' = s·cosφ − r·sinφ, since s = t × r and t × s = −r).
1981        let (ri, si) = match &twist_fractions {
1982            Some(fractions) => {
1983                let phi = twist_angle * fractions[station];
1984                let (sin_phi, cos_phi) = phi.sin_cos();
1985                let r = r_axes[station];
1986                let s = s_axes[station];
1987                (
1988                    r.scale(cos_phi).add(s.scale(sin_phi)),
1989                    s.scale(cos_phi).sub(r.scale(sin_phi)),
1990                )
1991            }
1992            None => (r_axes[station], s_axes[station]),
1993        };
1994        placed_axes.push((ri, si));
1995        let tk = tangents[station];
1996        let place = |x: Vec3| -> Vec3 {
1997            if rigid {
1998                let local = x.sub(path_start);
1999                station_origin
2000                    .add(tk.scale(local.dot(t0)))
2001                    .add(ri.scale(local.dot(r0)))
2002                    .add(si.scale(local.dot(s0)))
2003            } else {
2004                let local = x.sub(origin);
2005                station_origin
2006                    .add(ri.scale(local.dot(pu)))
2007                    .add(si.scale(local.dot(pv)))
2008            }
2009        };
2010        if rigid {
2011            // The section's own normal takes the same rotation (a direction, so
2012            // `R_k` without the translation), and the anchor origin's track is
2013            // the profile's bulk motion for the near-parallel guard.
2014            let n0 = placement.normal;
2015            section_normals.push(
2016                tk.scale(n0.dot(t0))
2017                    .add(ri.scale(n0.dot(r0)))
2018                    .add(si.scale(n0.dot(s0))),
2019            );
2020            anchor_track.push(place(origin));
2021            station_points.push(Vec::with_capacity(
2022                profile.iter().map(|c| c.control_points.len()).sum(),
2023            ));
2024        }
2025        let mut section = Vec::with_capacity(profile.len());
2026        for curve in profile {
2027            let control_points = curve
2028                .control_points
2029                .iter()
2030                .map(|point| {
2031                    let weight = point.w;
2032                    let euclidean = Vec3::new(point.x / weight, point.y / weight, point.z / weight);
2033                    let world = place(euclidean);
2034                    if rigid {
2035                        station_points[station].push(world);
2036                    }
2037                    Vec4 {
2038                        x: world.x * weight,
2039                        y: world.y * weight,
2040                        z: world.z * weight,
2041                        w: weight,
2042                    }
2043                })
2044                .collect();
2045            section.push(NurbsCurve::new(
2046                curve.degree,
2047                curve.knots.clone(),
2048                control_points,
2049            )?);
2050        }
2051        sections.push(section);
2052    }
2053
2054    // --- 4b. RIGID only: refuse a sweep whose sections fold through each other.
2055    //
2056    // TRANSPLANT keeps every section square to the path, so it can only fold
2057    // where the path's curvature radius drops below the profile's extent — the
2058    // case this builder has always documented as the caller's to avoid. RIGID
2059    // can fold for a second, much more reachable reason: a turning path carries
2060    // the profile around the TURN'S OWN AXIS, so a profile straddling that axis
2061    // has one half advancing while the other retreats. That is not a tolerance
2062    // question, it is a sign question, so it is answered exactly.
2063    //
2064    // The placement map is AFFINE in `x`, so the advance is affine in `x` too and
2065    // its extremes over the section lie on the control points. Checking the hull
2066    // therefore bounds the advance over the curves AND over the region they
2067    // enclose — this is a proof, not a sampling.
2068    //
2069    // ON A RING the guard walks one step FURTHER: the wrap from the last station
2070    // back to station 0 is a step the sections take like any other, and it is the
2071    // one step an open run does not have. A ring whose wrap step folded while
2072    // every interior step advanced would otherwise pass. On a ring whose twist
2073    // lands each curve `shift` along, the wrap step's control points arrive on
2074    // that curve's: the flat control-point list of station 0, turned left by the
2075    // points of the `shift` curves before it (curves `j` and `j + shift` carry the
2076    // same number, which is what makes the turn a relabelling).
2077    let step_count = if closure.is_ring() {
2078        stations
2079    } else {
2080        stations - 1
2081    };
2082    let next_station = |station: usize| (station + 1) % stations;
2083    if rigid {
2084        let landing_offset: usize = profile[..requested.shift]
2085            .iter()
2086            .map(|curve| curve.control_points.len())
2087            .sum();
2088        let mut most_positive = 0.0f64;
2089        let mut most_negative = 0.0f64;
2090        let mut longest_step = 0.0f64;
2091        for station in 0..step_count {
2092            let normal = section_normals[station];
2093            let arriving = &station_points[next_station(station)];
2094            let turn = if next_station(station) == 0 { landing_offset } else { 0 };
2095            for (before, after) in station_points[station]
2096                .iter()
2097                .zip(arriving[turn..].iter().chain(&arriving[..turn]))
2098            {
2099                let step = after.sub(*before);
2100                longest_step = longest_step.max(step.length());
2101                let advance = step.dot(normal);
2102                most_positive = most_positive.max(advance);
2103                most_negative = most_negative.min(advance);
2104            }
2105        }
2106        // Scale-relative, and far below any geometric tolerance: this band only
2107        // has to separate a sign from floating noise, never a small motion from
2108        // a large one.
2109        let noise = 1e-9 * longest_step.max(tolerance);
2110        if most_positive > noise && most_negative < -noise {
2111            return Err(
2112                "sweepSolid: the profile sweeps back through itself — part of it advances along \
2113                 the path while part of it retreats. A turning path carries the profile around \
2114                 the turn's own axis, so a profile that STRADDLES that axis folds into itself; \
2115                 move the profile clear of the axis, or sweep the run in pieces"
2116                    .into(),
2117            );
2118        }
2119        if most_positive <= noise && most_negative >= -noise {
2120            return Err(
2121                "sweepSolid: the path does not advance through the profile — it runs inside the \
2122                 profile plane, so the sweep encloses no volume"
2123                    .into(),
2124            );
2125        }
2126        // The near-parallel case, on the anchor origin's own track: the same 0.1
2127        // threshold `extrude_profile_brep` refuses a sliver at, so a STRAIGHT
2128        // path is accepted here exactly when the translational sweep of the same
2129        // profile and vector is. A station whose origin barely moves is skipped —
2130        // an origin sitting on the turn axis means the profile straddles it, and
2131        // the fold-back guard above has already spoken.
2132        for station in 0..step_count {
2133            let step = anchor_track[next_station(station)].sub(anchor_track[station]);
2134            let length = step.length();
2135            if length <= noise {
2136                continue;
2137            }
2138            if step.dot(section_normals[station]).abs() / length < 0.1 {
2139                return Err(
2140                    "sweepSolid: the path is nearly parallel to the profile plane, which sweeps a \
2141                     sliver rather than a solid"
2142                        .into(),
2143                );
2144            }
2145        }
2146    }
2147
2148    // --- 4c. TRANSPLANT: refuse a path that bends TIGHTER than the section
2149    //         reaches, before the loft rather than after it.
2150    //
2151    // A section square to the path, offset `d` towards the centre of the
2152    // osculating circle, advances `(1 − d·κ)` per unit of path: at `d·κ = 1` it
2153    // stands still and past it it RETREATS, so consecutive sections cross and
2154    // the wall the loft skins through them passes through itself. That is the
2155    // same `ρκ > 1` the fillet's canal surface folds at, and the same defect the
2156    // soundness acceptance catches at the loft's exit — but caught THERE it can
2157    // only name a face and two parameters, where here it can name the station,
2158    // the bend and the section that does not fit through it.
2159    //
2160    // RIGID is not measured here: it does not keep the section square to the
2161    // path, so `ρκ` is not its criterion, and §4b above is its guard.
2162    //
2163    // The reach is a WITNESS and the curvature is the PATH's own, so this refuses
2164    // only where a real section point really does turn back. What it does not see
2165    // is a bend that peaks BETWEEN two stations — the loft interpolates those,
2166    // and the acceptance at the exit remains the backstop.
2167    if !rigid && std::env::var("BREP_SWEEP_BEND_GUARD").as_deref() != Ok("0") {
2168        let offsets = section_offsets(profile, &placement)?;
2169        let profile_stations = bend_stations(points, curvatures, &placed_axes, &offsets);
2170        if let Some(refusal) = tight_bend_refusal(&profile_stations, stations) {
2171            // The ENVELOPE lane builds this shape exactly for ONE configuration
2172            // (a circle on one circular planar arc), and a caller that lands here
2173            // has missed it. `envelope_gap` is that lane's own measurement of
2174            // WHICH half was missed, carried down from the entry that has the
2175            // path to measure — so the refusal names the boundary of what is
2176            // built rather than only the bend that is refused.
2177            return Err(match envelope_gap {
2178                Some(gap) if !gap.is_empty() => format!(
2179                    "{refusal}. The swept envelope trimmed at its own self-intersection is built \
2180                     for a CIRCLE carried along one circular planar arc, and this sweep is not \
2181                     that: {gap}"
2182                ),
2183                _ => refusal,
2184            });
2185        }
2186    }
2187
2188    // --- 5. Loft through the swept sections.
2189    //
2190    //   OPEN: side walls plus two planar end caps.
2191    //
2192    //   RING: the CLOSED loft — the sections form a cycle, the last flowing back
2193    //   into the first through the same C² cyclic interpolation every interior
2194    //   span gets, so the seam is not a weld. One face per section curve and NO
2195    //   CAPS, which is the whole point: a ring has no ends to cap, and capping it
2196    //   would be two faces inside the material. The topology it emits is
2197    //   `c` corner vertices, `c` corner rings, `c` seam edges and `c` faces, so
2198    //   `V − E + F = c − 2c + c = 0` — the Euler characteristic of a torus, which
2199    //   is what a section carried once round a loop encloses.
2200    //   A RING whose requested twist lands each curve `shift` along skins each
2201    //   wall from its own curve into that one (`loft_profile_brep_closed_shifted`),
2202    //   so there is still one wall per section curve.
2203    match closure {
2204        RingFrame::Planar(_) | RingFrame::Spatial => {
2205            crate::loft_profile_brep_closed_shifted(&sections, requested.shift).map_err(|error| {
2206                format!("sweepSolid: closed loft through swept sections failed: {error}")
2207            })
2208        }
2209        RingFrame::Open => loft_profile_brep(&sections)
2210            .map_err(|error| format!("sweepSolid: loft through swept sections failed: {error}")),
2211    }
2212}
2213
2214/// Helical sweep (§5.7): sweep a CLOSED PLANAR profile loop along a helix of
2215/// `helix_radius` about the axis through `axis_origin` with direction
2216/// `axis_direction`, rising `pitch` per revolution for `turns` revolutions.
2217///
2218/// Sample a helix — `turns` revolutions about `axis_direction` from
2219/// `axis_origin`, starting at angle `start_angle` (radians, measured from the
2220/// axis frame's `u`) and rising `pitch` per turn, the radius running linearly
2221/// from `start_radius` to `end_radius` — uniformly in angle: 64 stations per
2222/// turn (at least 9, at most 4097 in total, so beyond 64 turns the per-turn
2223/// density thins). Returns the points and their chord parameters `s ∈ [0,1]`
2224/// (uniform-in-angle IS chord-length for a constant-radius helix, which is
2225/// what the averaged-knot interpolation assumes). `left_handed` winds the
2226/// angle the other way. The axis frame is `w = axis`, `u` = the component of
2227/// `reference` perpendicular to the axis (so angle zero points AT the
2228/// reference — the helix feature's local +X, or its start point), falling
2229/// back to `w.perpendicular()` when there is no usable reference, and
2230/// `v = w × u`. Shared by [`sweep_profile_helix`] and the HX feature, so the
2231/// coil a sweep builds and the edge a helix feature publishes agree exactly.
2232///
2233/// Errors on a degenerate axis, non-finite inputs, `turns ≤ 0`, more than 256
2234/// turns, a negative radius or pitch, or a helix that is a single point (both
2235/// radii AND the pitch zero).
2236#[allow(clippy::too_many_arguments)]
2237pub fn helix_sample_points(
2238    axis_origin: Vec3,
2239    axis_direction: Vec3,
2240    reference: Option<Vec3>,
2241    start_radius: f64,
2242    end_radius: f64,
2243    pitch: f64,
2244    turns: f64,
2245    start_angle: f64,
2246    left_handed: bool,
2247) -> Result<(Vec<Vec3>, Vec<f64>), String> {
2248    use std::f64::consts::TAU;
2249    let w = axis_direction
2250        .normalized()
2251        .map_err(|_| "helix: axis direction is degenerate".to_string())?;
2252    for (name, value) in [
2253        ("radius", start_radius),
2254        ("end radius", end_radius),
2255        ("pitch", pitch),
2256        ("turns", turns),
2257        ("start angle", start_angle),
2258    ] {
2259        if !value.is_finite() {
2260            return Err(format!("helix: {name} must be a finite number"));
2261        }
2262    }
2263    if start_radius < 0.0 || end_radius < 0.0 {
2264        return Err("helix: radius must not be negative".into());
2265    }
2266    if pitch < 0.0 {
2267        return Err("helix: pitch must not be negative".into());
2268    }
2269    if turns <= 0.0 {
2270        return Err("helix: turns must be positive".into());
2271    }
2272    if turns > 256.0 {
2273        return Err("helix: turns must be at most 256".into());
2274    }
2275    if start_radius == 0.0 && end_radius == 0.0 && pitch == 0.0 {
2276        return Err("helix: zero radius and zero pitch describe a single point".into());
2277    }
2278    let u = match reference {
2279        Some(reference) => {
2280            let radial = reference.sub(w.scale(reference.dot(w)));
2281            radial.normalized().or_else(|_| w.perpendicular())?
2282        }
2283        None => w.perpendicular()?,
2284    };
2285    let v = w.cross(u).normalized()?;
2286    let total_angle = turns * TAU;
2287    let height = pitch * turns;
2288    let sign = if left_handed { -1.0 } else { 1.0 };
2289    let count = ((turns * 64.0).ceil() as usize + 1).clamp(9, 4097);
2290    let mut points = Vec::with_capacity(count);
2291    let mut parameters = Vec::with_capacity(count);
2292    for index in 0..count {
2293        let s = index as f64 / (count - 1) as f64;
2294        let theta = start_angle + sign * total_angle * s;
2295        let radius = start_radius + (end_radius - start_radius) * s;
2296        points.push(
2297            axis_origin
2298                .add(u.scale(radius * theta.cos()))
2299                .add(v.scale(radius * theta.sin()))
2300                .add(w.scale(height * s)),
2301        );
2302        parameters.push(s);
2303    }
2304    Ok((points, parameters))
2305}
2306
2307/// [`helix_sample_points`] fitted as ONE cubic curve (global interpolation —
2308/// `interpolate_curve`, the machinery every fitted path here uses). Cubic
2309/// interpolation error on a circle of radius R with node spacing Δθ is
2310/// ≈ R·Δθ⁴/384: ~2·10⁻⁷·R at 64 stations per turn.
2311#[allow(clippy::too_many_arguments)]
2312pub fn fit_helix_curve(
2313    axis_origin: Vec3,
2314    axis_direction: Vec3,
2315    reference: Option<Vec3>,
2316    start_radius: f64,
2317    end_radius: f64,
2318    pitch: f64,
2319    turns: f64,
2320    start_angle: f64,
2321    left_handed: bool,
2322) -> Result<NurbsCurve, String> {
2323    let (points, parameters) = helix_sample_points(
2324        axis_origin,
2325        axis_direction,
2326        reference,
2327        start_radius,
2328        end_radius,
2329        pitch,
2330        turns,
2331        start_angle,
2332        left_handed,
2333    )?;
2334    interpolate_curve(&points, 3, &parameters)
2335}
2336
2337/// The helix is transcendental, not exactly NURBS-representable, so the path
2338/// is FIT: points (R·cosθ, R·sinθ, pitch·θ/2π) in the axis frame are sampled
2339/// uniformly in θ and globally interpolated with a cubic (`interpolate_curve`
2340/// — the same machinery every other fitted path here uses).  DENSITY: 64
2341/// samples per turn, capped at 1025 total nodes.  Cubic interpolation error
2342/// on a circle of radius R with node spacing Δθ is ≈ R·Δθ⁴/384: ~2·10⁻⁷·R at
2343/// 64/turn, and still ~6·10⁻⁵·R at the cap's worst case (16/turn at the
2344/// 64-turn limit) — orders below any profile a caller could sweep without
2345/// self-intersecting.  Uniform-in-θ parameters are exact chord-length for a
2346/// helix (constant speed), which is what the averaged-knot interpolation
2347/// assumes.
2348///
2349/// The fitted path then drives the EXISTING path-sweep core with 32 stations
2350/// per turn (min 32, capped at 1024 — the loft's dense interpolation solve is
2351/// O(stations³) per control column, so the cap trades per-turn density, never
2352/// correctness, at high turn counts).
2353///
2354/// GUARDS beyond the path sweep's own: the path sweep documents
2355/// self-intersection as caller responsibility, so the helix variant — which
2356/// knows its curvature analytically — rejects the two garbage modes itself:
2357///   • fold-over: the helix curvature radius (R² + c²)/R (c = pitch/2π) must
2358///     exceed the profile's max extent about its centroid, or the tube folds
2359///     through itself on the inner side (the torus tube-radius > major-radius
2360///     failure, pitch-relaxed);
2361///   • coil collision (turns ≥ 1): the normal gap between consecutive coils,
2362///     pitch·2πR/√((2πR)² + pitch²), must exceed the profile diameter.
2363/// Both use the profile's max sample distance from its centroid — conservative
2364/// for asymmetric profiles (extent in a harmless direction still counts), but
2365/// a false reject beats silent garbage.
2366///
2367/// SHALLOW COILS BUILD. They did not until 2026-09-08: the loft measured its
2368/// advance by the END-TO-END chord, which on a whole-turn helix is purely
2369/// AXIAL while both caps face nearly tangentially, so `|n·chord|` fell under
2370/// the loft's 0.1 sliver band and everything with pitch ≲ 0.63·R — most
2371/// springs — was refused. The loft measures the advance at each END now
2372/// (`loft_topology::basic::advance_from`), which is the direction the tube
2373/// actually grows through the cap it is orienting, so the band only ever
2374/// catches a real sliver.
2375pub fn sweep_profile_helix(
2376    profile: &[NurbsCurve],
2377    axis_origin: Vec3,
2378    axis_direction: Vec3,
2379    helix_radius: f64,
2380    pitch: f64,
2381    turns: f64,
2382    name: Option<&str>,
2383) -> Result<BrepSolid, String> {
2384    use std::f64::consts::TAU;
2385
2386    // --- 1. Validate the helix parameters with honest errors.
2387    let w = axis_direction
2388        .normalized()
2389        .map_err(|_| "sweep_profile_helix: axis direction is degenerate".to_string())?;
2390    if !(helix_radius.is_finite() && helix_radius > 0.0) {
2391        return Err("sweep_profile_helix: helix radius must be positive".into());
2392    }
2393    if !(pitch.is_finite() && pitch > 0.0) {
2394        return Err("sweep_profile_helix: pitch must be positive".into());
2395    }
2396    if !(turns.is_finite() && turns > 0.0) {
2397        return Err("sweep_profile_helix: turns must be positive".into());
2398    }
2399    // 64 turns bounds the loft's O(stations³) interpolation solve; beyond it
2400    // the station cap would silently degrade per-turn density anyway.
2401    if turns > 64.0 {
2402        return Err("sweep_profile_helix: turns must be at most 64".into());
2403    }
2404
2405    // --- 2. Profile extent about its centroid, sampled exactly like the path
2406    //        sweep derives its placement origin (boundary-sample mean), so the
2407    //        extent is measured about the point that actually rides the path.
2408    let mut samples = Vec::new();
2409    for curve in profile {
2410        let [start, end] = curve.domain()?;
2411        for sample in 0..16 {
2412            samples.push(curve.evaluate(start + (end - start) * sample as f64 / 16.0)?);
2413        }
2414    }
2415    if !samples.is_empty() {
2416        let mut centroid = Vec3::default();
2417        for point in &samples {
2418            centroid = centroid.add(*point);
2419        }
2420        let origin = centroid.scale(1.0 / samples.len() as f64);
2421        let extent = samples
2422            .iter()
2423            .map(|point| point.sub(origin).length())
2424            .fold(0.0, f64::max);
2425        let c = pitch / TAU; // axial rise per radian
2426                             // Fold-over: profile reaches past the helix's center of curvature.
2427        let curvature_radius = (helix_radius * helix_radius + c * c) / helix_radius;
2428        if extent >= curvature_radius {
2429            return Err(format!(
2430                "sweep_profile_helix: profile extent {extent:.6} reaches the helix \
2431                 curvature radius {curvature_radius:.6}; the tube would fold through \
2432                 itself — increase the helix radius or pitch, or shrink the profile"
2433            ));
2434        }
2435        // Coil collision: only possible once the sweep spans a full revolution.
2436        if turns >= 1.0 {
2437            let circumference = TAU * helix_radius;
2438            let turn_length = (circumference * circumference + pitch * pitch).sqrt();
2439            let coil_gap = pitch * circumference / turn_length;
2440            if coil_gap <= 2.0 * extent {
2441                return Err(format!(
2442                    "sweep_profile_helix: consecutive turns would self-intersect — \
2443                     coil gap {coil_gap:.6} does not clear the profile diameter {:.6}; \
2444                     increase the pitch or shrink the profile",
2445                    2.0 * extent
2446                ));
2447            }
2448        }
2449    }
2450
2451    // --- 3. Fit the helical path (see the density rationale in the fn docs) —
2452    //        the SAME sampler + fit the HX feature publishes as its edge.
2453    let path = fit_helix_curve(
2454        axis_origin,
2455        w,
2456        None,
2457        helix_radius,
2458        helix_radius,
2459        pitch,
2460        turns,
2461        0.0,
2462        false,
2463    )
2464    .map_err(|error| format!("sweep_profile_helix: helix path fit failed: {error}"))?;
2465
2466    // --- 4. Drive the existing sweep core; its (or the loft's) failures
2467    //        propagate with helix context prepended.
2468    let stations = ((turns * 32.0).ceil() as usize).clamp(32, 1024);
2469    // The helix builder places the profile ON the fitted helix (the coil's
2470    // section is drawn once and carried round), which is `Transplant` by
2471    // construction — the radius comes from the helix argument, not the profile's
2472    // own position.
2473    sweep_profile_along_path_stations(
2474        profile,
2475        &path,
2476        name,
2477        stations,
2478        0.0,
2479        None,
2480        SectionPlacement::Transplant,
2481    )
2482        .map_err(|error| format!("sweep_profile_helix: {error}"))
2483}