Skip to main content

brep_kernel/blending/blend/chain/
closed.rs

1use super::*;
2
3/// Piece of a support rim assigned to one mate face.
4pub(super) struct RimPiece {
5    pub(super) face_id: u64,
6    pub(super) loop_index: usize,
7    /// Global-parameter window of this piece.
8    pub(super) window: [f64; 2],
9    pub(super) edge_id: u64,
10}
11
12/// Blend a chain of conjugated edges with one rolling-ball blend face
13/// (Golovanov §6.9.5: all conjugated edges processed together).  A CLOSED
14/// chain (stadium rim, T-pipe saddle) welds into a periodic blend; an OPEN
15/// chain (line->arc->line capped by end faces) rolls a clamped blend that
16/// terminates in a transverse edge on each end face.
17pub fn blend_smooth_chain(
18    solid: &BrepSolid,
19    seed_edge_id: u64,
20    radius: f64,
21    chamfer: bool,
22    name: Option<&str>,
23) -> Result<BrepSolid, String> {
24    if !(radius > 0.0) || !radius.is_finite() {
25        return Err("blend: radius must be positive".into());
26    }
27    let chain = collect_smooth_chain(solid, seed_edge_id)?;
28    if chain.segments.len() < 2 {
29        return Err("blend: chain collapsed to a single segment".into());
30    }
31    if chain.closed {
32        blend_closed_smooth_chain(solid, &chain.segments, radius, chamfer, name)
33    } else {
34        blend_open_smooth_chain(solid, &chain, radius, chamfer, name)
35    }
36}
37
38/// Blend a CLOSED chain of conjugated edges with one rolling-ball blend
39/// face (Golovanov §6.9.5: all conjugated edges processed together).
40fn blend_closed_smooth_chain(
41    solid: &BrepSolid,
42    segments: &[ChainSegment<'_>],
43    radius: f64,
44    chamfer: bool,
45    name: Option<&str>,
46) -> Result<BrepSolid, String> {
47    let samples = march_chain(segments, radius, false)?;
48
49    // Global rows from the in-segment samples (wrapped-overlap closed fit).
50    let mut global: Vec<(f64, &ChainSample)> = samples
51        .iter()
52        .filter(|sample| (0..CHAIN_PER_SEGMENT as isize).contains(&sample.position))
53        .map(|sample| (sample.parameter.rem_euclid(1.0), sample))
54        .collect();
55    global.sort_by(|a, b| a.0.total_cmp(&b.0));
56    let count = global.len();
57    let degree = FIT_DEGREE;
58    let mut extended_params = Vec::new();
59    let mut samples_cr = Vec::new();
60    let mut samples_mid = Vec::new();
61    let mut samples_cs = Vec::new();
62    let mut push = |station: &Station, parameter: f64| {
63        extended_params.push(parameter);
64        samples_cr.push(Vec4::from_point(station.p1, 1.0));
65        samples_cs.push(Vec4::from_point(station.p2, 1.0));
66        samples_mid.push(Vec4 {
67            x: station.apex.x * station.weight,
68            y: station.apex.y * station.weight,
69            z: station.apex.z * station.weight,
70            w: station.weight,
71        });
72    };
73    for offset in (1..=degree).rev() {
74        let (parameter, sample) = &global[count - offset];
75        push(&sample.station, parameter - 1.0);
76    }
77    for (parameter, sample) in &global {
78        push(&sample.station, *parameter);
79    }
80    // Close the loop: repeat the first station at parameter 1, then the
81    // wrap continuation.
82    push(&global[0].1.station, global[0].0 + 1.0);
83    for offset in 1..=degree {
84        let (parameter, sample) = &global[offset];
85        push(&sample.station, parameter + 1.0);
86    }
87    let low = extended_params[0];
88    let high = *extended_params.last().unwrap();
89    let range = high - low;
90    let normalized: Vec<f64> = extended_params
91        .iter()
92        .map(|parameter| (parameter - low) / range)
93        .collect();
94    let seam_low = (0.0 - low) / range;
95    let seam_high = (1.0 - low) / range;
96    let fit_row = |row: &[Vec4]| -> Result<NurbsCurve, String> {
97        let curve = fit::interpolate_homogeneous(row, degree, &normalized)?;
98        let (_, tail) = curve.split(seam_low)?;
99        let (middle, _) = tail.split(seam_high)?;
100        Ok(middle)
101    };
102    let cr = fit_row(&samples_cr)?;
103    let cs = fit_row(&samples_cs)?;
104    let mid = if chamfer {
105        None
106    } else {
107        Some(fit_row(&samples_mid)?)
108    };
109    let u_domain = cr.domain()?;
110    let rows_u = cr.control_points.len();
111    let mut control = Vec::with_capacity(rows_u);
112    for index in 0..rows_u {
113        let mut column = vec![cr.control_points[index]];
114        if let Some(mid) = &mid {
115            column.push(mid.control_points[index]);
116        }
117        column.push(cs.control_points[index]);
118        control.push(column);
119    }
120    let last = rows_u - 1;
121    for column_index in 0..control[0].len() {
122        control[last][column_index] = control[0][column_index];
123    }
124    let (degree_v, knots_v) = if chamfer {
125        (1usize, vec![0.0, 0.0, 1.0, 1.0])
126    } else {
127        (2usize, vec![0.0, 0.0, 0.0, 1.0, 1.0, 1.0])
128    };
129    let surface = NurbsSurface::new(degree, degree_v, cr.knots.clone(), knots_v, control)?;
130
131    // Single-face sides (e.g. one cap around the whole rim) get ONE uv
132    // fit in the rows' own normalized space — identical parameterization
133    // to cr/cs, so the closed support piece and its pcurve agree exactly.
134    let single_face = |side: usize| -> bool {
135        let first_id = if side == 0 {
136            segments[0].first.face.id
137        } else {
138            segments[0].second.face.id
139        };
140        segments.iter().all(|segment| {
141            let id = if side == 0 {
142                segment.first.face.id
143            } else {
144                segment.second.face.id
145            };
146            id == first_id
147        })
148    };
149    let mut whole_pcurves: [Option<NurbsCurve>; 2] = [None, None];
150    for side in 0..2 {
151        if !single_face(side) {
152            continue;
153        }
154        let mut row = Vec::new();
155        let select = |station: &Station| -> [f64; 2] {
156            if side == 0 {
157                station.uv1
158            } else {
159                station.uv2
160            }
161        };
162        let mut push_uv = |station: &Station| {
163            let uv = select(station);
164            row.push(Vec4::from_point(Vec3::new(uv[0], uv[1], 0.0), 1.0));
165        };
166        for offset in (1..=degree).rev() {
167            push_uv(&global[count - offset].1.station);
168        }
169        for (_, sample) in &global {
170            push_uv(&sample.station);
171        }
172        push_uv(&global[0].1.station);
173        for offset in 1..=degree {
174            push_uv(&global[offset].1.station);
175        }
176        whole_pcurves[side] = Some(fit_row(&row)?);
177    }
178
179    // Per-face pcurve fits over each segment's full sample window
180    // (in-segment + overshoot), in global parameters.
181    let mut pcurves1 = Vec::with_capacity(segments.len());
182    let mut pcurves2 = Vec::with_capacity(segments.len());
183    for segment in 0..segments.len() {
184        let mut window: Vec<&ChainSample> = samples
185            .iter()
186            .filter(|sample| sample.segment == segment)
187            .collect();
188        window.sort_by(|a, b| a.parameter.total_cmp(&b.parameter));
189        let params: Vec<f64> = window.iter().map(|sample| sample.parameter).collect();
190        let low = params[0];
191        let high = *params.last().unwrap();
192        let normalized: Vec<f64> = params
193            .iter()
194            .map(|parameter| (parameter - low) / (high - low))
195            .collect();
196        let fit_uv = |select: &dyn Fn(&Station) -> [f64; 2]| -> Result<NurbsCurve, String> {
197            let points: Vec<Vec4> = window
198                .iter()
199                .map(|sample| {
200                    let uv = select(&sample.station);
201                    Vec4::from_point(Vec3::new(uv[0], uv[1], 0.0), 1.0)
202                })
203                .collect();
204            let curve = fit::interpolate_homogeneous(&points, degree, &normalized)?;
205            // Re-express in GLOBAL parameters: the curve's [0,1] domain
206            // corresponds to [low, high]; keep as-is and remember the
207            // affine map through the window bounds.
208            Ok(curve)
209        };
210        pcurves1.push((low, high, fit_uv(&|station| station.uv1)?));
211        pcurves2.push((low, high, fit_uv(&|station| station.uv2)?));
212    }
213
214    chain_surgery(
215        solid,
216        segments,
217        ChainRows {
218            surface,
219            cr,
220            cs,
221            u_domain,
222            fit_low: low,
223            fit_range: range,
224            pcurves1,
225            pcurves2,
226            whole_pcurves,
227        },
228        name,
229    )
230}
231
232pub(super) struct ChainRows {
233    pub(super) surface: NurbsSurface,
234    pub(super) cr: NurbsCurve,
235    pub(super) cs: NurbsCurve,
236    pub(super) u_domain: [f64; 2],
237    /// Affine map from the rows' fit space to GLOBAL chain parameters:
238    /// global = fit_low + p * fit_range (the global fit normalised its
239    /// wrap-extended parameters onto [0, 1] before interpolating).
240    pub(super) fit_low: f64,
241    pub(super) fit_range: f64,
242    /// Per-segment (window_low, window_high, uv curve over [0,1]) in
243    /// GLOBAL parameters.
244    pub(super) pcurves1: Vec<(f64, f64, NurbsCurve)>,
245    pub(super) pcurves2: Vec<(f64, f64, NurbsCurve)>,
246    /// Whole-chain uv fits (rows' fit space) for single-face sides.
247    pub(super) whole_pcurves: [Option<NurbsCurve>; 2],
248}
249
250impl ChainRows {
251    pub(super) fn to_global(&self, fit_parameter: f64) -> f64 {
252        self.fit_low + fit_parameter * self.fit_range
253    }
254}
255
256/// Extract the pcurve portion for a global window [a, b] from a
257/// segment's fitted uv curve, re-parameterised so its domain maps
258/// affinely onto the portion (matching the split support edge).
259pub(super) fn pcurve_portion(fitted: &(f64, f64, NurbsCurve), a: f64, b: f64) -> Result<NurbsCurve, String> {
260    let (low, high, curve) = fitted;
261    // The chain parameterisation wraps with period 1; pick the period
262    // copy of the window that overlaps this segment's fit range.
263    let to_local = |value: f64| {
264        let mut best = value;
265        for candidate in [value, value + 1.0, value - 1.0] {
266            if candidate >= low - 1e-9 && candidate <= high + 1e-9 {
267                best = candidate;
268                break;
269            }
270        }
271        ((best - low) / (high - low)).clamp(0.0, 1.0)
272    };
273    let mut local_a = to_local(a);
274    let mut local_b = to_local(b);
275    if local_a > local_b {
276        std::mem::swap(&mut local_a, &mut local_b);
277    }
278    let epsilon = 1e-9;
279    let (_, tail) = if local_a > epsilon {
280        curve.split(local_a)?
281    } else {
282        (curve.clone(), curve.clone())
283    };
284    let tail = if local_a > epsilon {
285        tail
286    } else {
287        curve.clone()
288    };
289    let portion = if local_b < 1.0 - epsilon {
290        tail.split(local_b)?.0
291    } else {
292        tail
293    };
294    Ok(portion)
295}
296
297/// The single boundary edge of a side's mate face(s) incident to a chain
298/// vertex, EXCLUDING the two chain edges meeting there — i.e. the SEAM
299/// (same face on both segments) or the SPOKE (face changes) that the
300/// support curve must cross.  `None` when the chain simply flows through
301/// (same face, no seam — e.g. a stadium cap rim vertex).
302pub(super) fn cross_edge_at(
303    solid: &BrepSolid,
304    face_before: &FaceRecord,
305    loop_before: usize,
306    face_after: &FaceRecord,
307    loop_after: usize,
308    vertex: u64,
309    before_edge: u64,
310    after_edge: u64,
311) -> Result<Option<u64>, String> {
312    let mut found: Option<u64> = None;
313    let mut consider = |face: &FaceRecord, loop_index: usize| -> Result<(), String> {
314        for coedge in &face.loops[loop_index].coedges {
315            if coedge.edge_id == before_edge || coedge.edge_id == after_edge {
316                continue;
317            }
318            let edge = solid
319                .edges
320                .iter()
321                .find(|edge| edge.id == coedge.edge_id)
322                .ok_or("blend: loop references missing edge")?;
323            if edge.start_vertex_id == vertex || edge.end_vertex_id == vertex {
324                if found.is_some() && found != Some(edge.id) {
325                    return Err("blend: multiple cross edges at a chain vertex".into());
326                }
327                found = Some(edge.id);
328            }
329        }
330        Ok(())
331    };
332    consider(face_before, loop_before)?;
333    if !(face_after.id == face_before.id && loop_after == loop_before) {
334        consider(face_after, loop_after)?;
335    }
336    Ok(found)
337}
338
339/// The mate-face pcurve of a support PIECE, built by projecting the piece's
340/// 3D curve onto the (analytic) mate surface and fitting.  Because a chain
341/// support rides EXACTLY on its mate carrier, the projection is exact; the
342/// u track is unwrapped so it never jumps a period across the carrier seam,
343/// and crossing endpoints are pinned to the seam meridian (u0 or u1) so the
344/// piece and the trimmed seam edge close the loop.
345pub(super) fn project_piece_pcurve(
346    piece: &NurbsCurve,
347    surface: &NurbsSurface,
348    start_is_crossing: bool,
349    end_is_crossing: bool,
350) -> Result<NurbsCurve, String> {
351    let [u0, u1] = surface.domain_u()?;
352    let period = u1 - u0;
353    let [d0, d1] = piece.domain()?;
354    // Dense reference projection; the pcurve is fit from an adaptively
355    // thinned subset kept as coarse as tolerance allows.
356    const DENSE: usize = 96;
357    let mut proj_u = Vec::with_capacity(DENSE + 1);
358    let mut proj_v = Vec::with_capacity(DENSE + 1);
359    let mut point3 = Vec::with_capacity(DENSE + 1);
360    let mut params = Vec::with_capacity(DENSE + 1);
361    let mut previous_u: Option<f64> = None;
362    for section in 0..=DENSE {
363        let t = d0 + (d1 - d0) * section as f64 / DENSE as f64;
364        let point = piece.evaluate(t)?;
365        let projection = crate::project_point_to_surface(surface, point)?;
366        let mut u = projection.u;
367        if let Some(previous) = previous_u {
368            while u - previous > period * 0.5 {
369                u -= period;
370            }
371            while previous - u > period * 0.5 {
372                u += period;
373            }
374        }
375        previous_u = Some(u);
376        proj_u.push(u);
377        proj_v.push(projection.v);
378        point3.push(point);
379        params.push(section as f64 / DENSE as f64);
380    }
381    // Interior mean decides which meridian a crossing endpoint sits on, and
382    // whether the whole track needs a full-period shift into the domain.
383    let mean_u: f64 = proj_u[1..DENSE].iter().sum::<f64>() / (DENSE - 1) as f64;
384    let seam_u = if mean_u - u0 > u1 - mean_u { u1 } else { u0 };
385    if start_is_crossing {
386        proj_u[0] = seam_u;
387    }
388    if end_is_crossing {
389        proj_u[DENSE] = seam_u;
390    }
391    let mut shift = 0.0;
392    while mean_u + shift > u1 + 1e-9 {
393        shift -= period;
394    }
395    while mean_u + shift < u0 - 1e-9 {
396        shift += period;
397    }
398    for u in proj_u.iter_mut() {
399        *u += shift;
400    }
401    // Coarsest interpolation whose fitted pcurve stays inside HALF the
402    // validator's pcurve/edge tolerance everywhere along the dense track.
403    let tolerance = 0.002;
404    let mut best: Option<NurbsCurve> = None;
405    for sections in [8usize, 12, 16, 24, 32, 48, 64, 96] {
406        if sections > DENSE {
407            break;
408        }
409        let indices: Vec<usize> = (0..=sections)
410            .map(|k| (k * DENSE / sections).min(DENSE))
411            .collect();
412        let points: Vec<Vec4> = indices
413            .iter()
414            .map(|&i| Vec4::from_point(Vec3::new(proj_u[i], proj_v[i], 0.0), 1.0))
415            .collect();
416        let knot_params: Vec<f64> = indices.iter().map(|&i| params[i]).collect();
417        let curve = fit::interpolate_homogeneous(&points, FIT_DEGREE, &knot_params)?;
418        let mut max_deviation: f64 = 0.0;
419        for i in 0..=DENSE {
420            let uv = curve.evaluate(params[i])?;
421            let on_surface = surface.evaluate(uv.x, uv.y)?;
422            max_deviation = max_deviation.max(on_surface.sub(point3[i]).length());
423        }
424        best = Some(curve);
425        if max_deviation <= tolerance {
426            break;
427        }
428    }
429    best.ok_or_else(|| "blend: piece pcurve projection failed".to_string())
430}