use crate::{KernelRefusal, KernelStage, OrRefuse};
use super::*;
fn rebuild_pcurve_preserving_branch(
surface: &NurbsSurface,
curve: &NurbsCurve,
previous: &NurbsCurve,
emitted: bool,
) -> Result<NurbsCurve, KernelRefusal> {
if emitted {
let [t0, t1] = curve.domain().or_refuse(KernelStage::Intersect, "domain")?;
let start = curve.evaluate(t0).or_refuse(KernelStage::Intersect, "evaluate")?;
let end = curve.evaluate(t1).or_refuse(KernelStage::Intersect, "evaluate")?;
return super::hermite::piece_pcurve_on_chart(previous, surface, start, end).or_refuse(KernelStage::Intersect, "piece_pcurve_on_chart");
}
let mut rebuilt = build_pcurve_on_surface_marched(surface, curve).or_refuse(KernelStage::Intersect, "build_pcurve_on_surface_marched")?;
if std::env::var("BREP_PCURVE_BRANCH_ALIGN").as_deref() == Ok("0") {
return Ok(rebuilt);
}
let (closed_u, closed_v) = surface.closed_directions().or_refuse(KernelStage::Intersect, "closed_directions")?;
if !closed_u && !closed_v {
return Ok(rebuilt);
}
let mean = |pcurve: &NurbsCurve| -> Option<(f64, f64)> {
let mut sum_u = 0.0;
let mut sum_v = 0.0;
let mut count = 0.0;
for point in &pcurve.control_points {
if point.w.abs() > 1e-300 {
sum_u += point.x / point.w;
sum_v += point.y / point.w;
count += 1.0;
}
}
(count > 0.0).then(|| (sum_u / count, sum_v / count))
};
let (Some((old_u, old_v)), Some((new_u, new_v))) = (mean(previous), mean(&rebuilt)) else {
return Ok(rebuilt);
};
let [u0, u1] = surface.domain_u().or_refuse(KernelStage::Intersect, "domain_u")?;
let [v0, v1] = surface.domain_v().or_refuse(KernelStage::Intersect, "domain_v")?;
let shift_u = if closed_u && u1 > u0 {
((old_u - new_u) / (u1 - u0)).round()
} else {
0.0
};
let shift_v = if closed_v && v1 > v0 {
((old_v - new_v) / (v1 - v0)).round()
} else {
0.0
};
if shift_u == 0.0 && shift_v == 0.0 {
return Ok(rebuilt);
}
for point in &mut rebuilt.control_points {
point.x += shift_u * (u1 - u0) * point.w;
point.y += shift_v * (v1 - v0) * point.w;
}
Ok(rebuilt)
}
#[derive(Clone, Copy, Debug, PartialEq, Eq, PartialOrd, Ord)]
pub(super) enum JunctionSource {
OperandVertex,
EdgePoint,
FacePoint,
}
#[derive(Default)]
pub(super) struct JunctionEvidence {
pub(super) radii: HashMap<u64, f64>,
pub(super) sources: HashMap<u64, JunctionSource>,
}
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
pub(super) enum JunctionCombineMode {
Off,
Check,
On,
}
pub(super) fn junction_combine_mode(options: &ImprintOptions) -> JunctionCombineMode {
match mode_switch!(live "BREP_JUNCTION_COMBINE", options.two_radius_junctions) {
crate::env_read::Mode::Off => JunctionCombineMode::Off,
crate::env_read::Mode::Check => JunctionCombineMode::Check,
crate::env_read::Mode::On => JunctionCombineMode::On,
}
}
const STRICT_FACTOR: f64 = 5.0;
fn junction_census_line(line: impl FnOnce() -> String) {
if let Ok(path) = std::env::var("BREP_JUNCTION_CENSUS") {
use std::io::Write;
if let Ok(mut file) = std::fs::OpenOptions::new().create(true).append(true).open(path) {
let label = std::env::var("BREP_CENSUS_LABEL").unwrap_or_default();
let _ = writeln!(file, "junction census: label {label} {}", line());
}
}
}
fn junction_incidence(
result: &ImprintResultRecord,
) -> Result<(HashMap<u64, HashSet<FaceKey>>, HashMap<u64, Vec<Vec3>>), KernelRefusal> {
let mut incident_faces = HashMap::<u64, HashSet<FaceKey>>::default();
let mut incident_points = HashMap::<u64, Vec<Vec3>>::default();
for piece in &result.pieces {
let [start, end] = piece.curve.domain().or_refuse(KernelStage::Intersect, "domain")?;
for (vertex_id, point) in [
(piece.start_vertex_id, piece.curve.evaluate(start).or_refuse(KernelStage::Intersect, "evaluate")?),
(piece.end_vertex_id, piece.curve.evaluate(end).or_refuse(KernelStage::Intersect, "evaluate")?),
] {
let entry = incident_faces.entry(vertex_id).or_default();
entry.extend(piece.support_faces);
incident_points.entry(vertex_id).or_default().push(point);
}
}
Ok((incident_faces, incident_points))
}
#[allow(clippy::too_many_arguments)]
fn two_radius_combine(
result: &mut ImprintResultRecord,
solid_a: &BrepSolid,
solid_b: &BrepSolid,
face_map: &HashMap<FaceKey, &FaceRecord>,
charts: &FaceCharts<'_>,
incident_faces: &HashMap<u64, HashSet<FaceKey>>,
tolerance: f64,
raw_extent: f64,
evidence: &JunctionEvidence,
emitted_pcurves: &HashSet<(u64, FaceKey)>,
mode: JunctionCombineMode,
anchored: &mut HashSet<u64>,
) -> Result<bool, KernelRefusal> {
let weld = assembler_weld(tolerance);
let (_, capture) = residual_merge_bands(raw_extent, tolerance);
let mut face_vertices: HashMap<FaceKey, Vec<(u64, Vec3)>> = HashMap::default();
for (operand, solid) in [(0u8, solid_a), (1u8, solid_b)] {
let points: HashMap<u64, Vec3> = solid.vertices.iter().map(|vertex| (vertex.id, vertex.point)).collect();
let ends: HashMap<u64, (u64, u64)> =
solid.edges.iter().map(|edge| (edge.id, (edge.start_vertex_id, edge.end_vertex_id))).collect();
for face in solid.shells.iter().flat_map(|shell| &shell.faces) {
let mut ids: Vec<u64> = Vec::new();
for coedge in face.loops.iter().flat_map(|loop_record| &loop_record.coedges) {
if let Some(&(start, end)) = ends.get(&coedge.edge_id) {
for id in [start, end] {
if !ids.contains(&id) {
ids.push(id);
}
}
}
}
let entries = ids.into_iter().filter_map(|id| points.get(&id).map(|point| (id, *point))).collect();
face_vertices.insert(FaceKey { operand, face_id: face.id }, entries);
}
}
let radius_of = |id: u64| evidence.radii.get(&id).copied().unwrap_or(0.0).max(weld);
let source_of = |id: u64| evidence.sources.get(&id).copied().unwrap_or(JunctionSource::FacePoint);
let piece_connected: HashSet<(u64, u64)> = result
.pieces
.iter()
.map(|piece| (piece.start_vertex_id.min(piece.end_vertex_id), piece.start_vertex_id.max(piece.end_vertex_id)))
.collect();
let verts: Vec<(u64, Vec3)> = result.vertices.iter().map(|vertex| (vertex.id, vertex.point)).collect();
let empty = HashSet::default();
let mut anchor: HashMap<u64, (u8, u64, Vec3)> = HashMap::default();
let mut nearest_junction: HashMap<u64, (f64, u64, f64)> = HashMap::default();
for &(id, point) in &verts {
let keys = incident_faces.get(&id).unwrap_or(&empty);
let radius = radius_of(id);
let mut nearest_vertex: Option<(f64, u8, u64, Vec3)> = None;
for key in keys {
for &(vertex_id, vertex_point) in face_vertices.get(key).map(Vec::as_slice).unwrap_or(&[]) {
let distance = vertex_point.sub(point).length();
if distance <= capture && nearest_vertex.is_none_or(|best| distance < best.0) {
nearest_vertex = Some((distance, key.operand, vertex_id, vertex_point));
}
}
}
let mut nearest: Option<(f64, u64, f64)> = None;
for &(other_id, other_point) in &verts {
if other_id == id || piece_connected.contains(&(id.min(other_id), id.max(other_id))) {
continue;
}
let other_keys = incident_faces.get(&other_id).unwrap_or(&empty);
if keys.is_disjoint(other_keys) {
continue;
}
let distance = other_point.sub(point).length();
if distance <= capture && nearest.is_none_or(|best| distance < best.0) {
nearest = Some((distance, other_id, STRICT_FACTOR * (radius + radius_of(other_id))));
}
}
let strict_vertex = STRICT_FACTOR * (radius + weld);
junction_census_line(|| {
format!(
"vertex={id} source={:?} radius={radius:.3e} faces={} capture={capture:.3e} operand_vertex={} junction={}",
source_of(id),
keys.len(),
match nearest_vertex {
Some((distance, operand, vertex_id, _)) => format!(
"{distance:.3e}/{operand}:{vertex_id} strict={strict_vertex:.3e} inside={}",
distance <= strict_vertex
),
None => "none".to_string(),
},
match nearest {
Some((distance, other_id, strict)) => format!(
"{distance:.3e}/{other_id}:{:?} strict={strict:.3e} inside={}",
source_of(other_id),
distance <= strict
),
None => "none".to_string(),
}
)
});
if let Some((distance, operand, vertex_id, vertex_point)) = nearest_vertex {
if distance > tolerance && distance <= strict_vertex {
anchor.insert(id, (operand, vertex_id, vertex_point));
}
}
if let Some(found) = nearest {
nearest_junction.insert(id, found);
}
}
if mode != JunctionCombineMode::On {
return Ok(false);
}
let mut parent: HashMap<u64, u64> = verts.iter().map(|(id, _)| (*id, *id)).collect();
let mut merged_point: HashMap<u64, Vec3> = HashMap::default();
let anchors = anchor.clone();
let effective = |id: u64| -> (JunctionSource, f64, Vec3) {
match anchors.get(&id) {
Some((_, _, point)) => (JunctionSource::OperandVertex, weld, *point),
None => (
source_of(id),
radius_of(id),
verts.iter().find(|(other, _)| *other == id).map(|(_, point)| *point).unwrap_or_default(),
),
}
};
let mut changed = false;
for &(id, _) in &verts {
let Some(&(distance, other_id, strict)) = nearest_junction.get(&id) else {
continue;
};
if distance <= tolerance || distance > strict {
continue;
}
let (ra, rb) = (junction_find(&mut parent, id), junction_find(&mut parent, other_id));
if ra == rb {
continue;
}
let (source_a, radius_a, point_a) = effective(ra);
let (source_b, radius_b, point_b) = effective(rb);
let keep_a = (source_a, radius_a, ra) <= (source_b, radius_b, rb);
let (rep, absorbed, point) = if keep_a { (ra, rb, point_a) } else { (rb, ra, point_b) };
parent.insert(absorbed, rep);
merged_point.insert(rep, point);
if anchor.contains_key(&absorbed) && !anchor.contains_key(&rep) {
if let Some(found) = anchor.get(&absorbed).copied() {
anchor.insert(rep, found);
}
}
if std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
eprintln!(
"[two-radius] v{rep} <- v{absorbed} sep={distance:.3e} strict={strict:.3e} sources {source_a:?}/{source_b:?} to=({:.7},{:.7},{:.7})",
point.x, point.y, point.z
);
}
changed = true;
}
let edge_curves: HashMap<(u8, u64), &EdgeRecord> = [(0u8, solid_a), (1u8, solid_b)]
.into_iter()
.flat_map(|(operand, solid)| solid.edges.iter().map(move |edge| ((operand, edge.id), edge)))
.collect();
for (&id, &(operand, vertex_id, point)) in &anchor {
let rep = junction_find(&mut parent, id);
merged_point.insert(rep, point);
anchored.insert(rep);
let strict = STRICT_FACTOR * (radius_of(id) + weld);
for record in result.edge_splits.iter_mut().filter(|record| record.operand == operand) {
let Some(edge) = edge_curves.get(&(operand, record.edge_id)) else {
continue;
};
if edge.start_vertex_id != vertex_id && edge.end_vertex_id != vertex_id {
continue;
}
record.parameters.retain(|¶meter| {
edge.curve
.evaluate(parameter)
.map(|split_point| split_point.sub(point).length() > strict)
.unwrap_or(true)
});
}
if std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
eprintln!(
"[two-radius] v{id} anchored to operand {operand} vertex {vertex_id} at ({:.7},{:.7},{:.7})",
point.x, point.y, point.z
);
}
changed = true;
}
if !changed {
return Ok(false);
}
let (_, sep_cap) = residual_merge_bands(raw_extent, tolerance);
relocate_junctions(result, &mut parent, &merged_point, face_map, charts, emitted_pcurves, (4.0 * sep_cap).max(1e-3), true)?;
Ok(true)
}
#[allow(clippy::too_many_arguments)]
pub(super) fn canonicalize_imprint_junctions(
result: &mut ImprintResultRecord,
solid_a: &BrepSolid,
solid_b: &BrepSolid,
charts: &FaceCharts<'_>,
tolerance: f64,
emitted_pcurves: &HashSet<(u64, FaceKey)>,
evidence: &JunctionEvidence,
combine: JunctionCombineMode,
residual_merge: bool,
) -> Result<(), KernelRefusal> {
let face_map = faces(solid_a, 0)
.into_iter()
.chain(faces(solid_b, 1))
.map(|face| (face.key(), face.face))
.collect::<HashMap<_, _>>();
let (mut incident_faces, mut incident_points) = junction_incidence(result)?;
let edge_curves: HashMap<(u8, u64), &EdgeRecord> = [(0u8, solid_a), (1u8, solid_b)]
.into_iter()
.flat_map(|(operand, solid)| solid.edges.iter().map(move |edge| ((operand, edge.id), edge)))
.collect();
let operand_vertices: Vec<(u8, Vec3)> = solid_a
.vertices
.iter()
.map(|vertex| (0u8, vertex.point))
.chain(solid_b.vertices.iter().map(|vertex| (1u8, vertex.point)))
.collect();
let raw_extent = raw_solid_extent(solid_a).max(raw_solid_extent(solid_b));
if residual_merge
&& std::env::var("BREP_IMPRINT_RESIDUAL_MERGE").as_deref() != Ok("0")
&& merge_residual_coincident_vertices(
result,
&face_map,
&edge_curves,
&operand_vertices,
charts,
&incident_faces,
tolerance,
raw_extent,
evidence,
emitted_pcurves,
)?
{
(incident_faces, incident_points) = junction_incidence(result)?;
}
let mut anchored: HashSet<u64> = HashSet::default();
if combine != JunctionCombineMode::Off
&& two_radius_combine(
result,
solid_a,
solid_b,
&face_map,
charts,
&incident_faces,
tolerance,
raw_extent,
evidence,
emitted_pcurves,
combine,
&mut anchored,
)?
{
(incident_faces, incident_points) = junction_incidence(result)?;
}
if std::env::var("BREP_RESIDUAL_MERGE_EDGE_PLACEMENT").as_deref() != Ok("0")
&& consolidate_fuzzy_piece_ends(result, &face_map, &edge_curves, charts, tolerance, evidence, emitted_pcurves)?
{
(incident_faces, incident_points) = junction_incidence(result)?;
}
let reference = junction_reference(result.vertices.iter().map(|vertex| vertex.point));
if std::env::var("BREP_DEBUG_CANONICAL_JUNCTIONS").is_ok() {
eprintln!("JUNCTION-REFERENCE point={reference:?} vertices={}", result.vertices.len());
}
let mut canonical_vertices = result
.vertices
.iter()
.filter(|vertex| {
let limit = tolerance * 10.0 * (1.0 + vertex.point.sub(reference).length());
incident_points.get(&vertex.id).is_some_and(|points| {
points
.iter()
.all(|point| point.sub(vertex.point).length() <= limit)
})
})
.map(|vertex| vertex.id)
.collect::<HashSet<_>>();
if std::env::var("BREP_DEBUG_CANONICAL_JUNCTIONS").is_ok() {
for vertex in &result.vertices {
let limit = tolerance * 10.0 * (1.0 + vertex.point.sub(reference).length());
let gap = incident_points.get(&vertex.id).map(|points| points.iter()
.map(|p| p.sub(vertex.point).length()).fold(0.0_f64, f64::max));
let old_limit = tolerance * 10.0 * (1.0 + vertex.point.length());
let old_member = gap.is_some_and(|gap| gap <= old_limit);
let member = canonical_vertices.contains(&vertex.id);
eprintln!("JUNCTION-MEMBERSHIP vertex={} point={:?} tolerance={tolerance:.17e} incident_gap={gap:?} limit={limit:.17e} canonical={member} old_limit={old_limit:.17e} old_canonical={old_member} flip={}", vertex.id, vertex.point, member != old_member);
if member != old_member {
eprintln!("JUNCTION-FLIP vertex={} faces={:?} pieces={:?}", vertex.id, incident_faces.get(&vertex.id), result.pieces.iter().filter(|piece| piece.start_vertex_id == vertex.id || piece.end_vertex_id == vertex.id).map(|piece| piece.id).collect::<Vec<_>>());
}
}
}
for vertex in &mut result.vertices {
if !canonical_vertices.contains(&vertex.id) {
continue;
}
if anchored.contains(&vertex.id) {
continue;
}
let Some(keys) = incident_faces.get(&vertex.id) else {
continue;
};
let surfaces = keys
.iter()
.filter_map(|key| face_map.get(key).map(|face| &face.surface))
.collect::<Vec<_>>();
if surfaces.len() < 2 {
continue;
}
let residual = |point: Vec3| -> Result<f64, KernelRefusal> {
surfaces
.iter()
.map(|surface| Ok(project_point_to_surface(surface, point).or_refuse(KernelStage::Intersect, "project_point_to_surface")?.distance))
.try_fold(0.0f64, |maximum, distance| {
distance.map(|distance| maximum.max(distance))
})
};
let mut current = vertex.point;
let mut best = current;
let mut best_residual = residual(current)?;
for _ in 0..60 {
let projected = surfaces
.iter()
.map(|surface| project_point_to_surface(surface, current).map(|value| value.point))
.collect::<Result<Vec<_>, _>>().or_refuse(KernelStage::Intersect, "csg.imprint.junctions")?;
let candidate = projected
.iter()
.copied()
.fold(Vec3::default(), Vec3::add)
.scale(1.0 / projected.len() as f64);
let candidate_residual = residual(candidate)?;
if candidate_residual < best_residual {
best = candidate;
best_residual = candidate_residual;
}
if candidate.sub(current).length() <= tolerance.max(1e-12) {
break;
}
current = candidate;
}
let maximum_move = tolerance * 100.0 * (1.0 + vertex.point.sub(reference).length());
if std::env::var("BREP_DEBUG_CANONICAL_JUNCTIONS").is_ok() {
eprintln!("JUNCTION-POLISH vertex={} point={:?} candidate={best:?} tolerance={tolerance:.17e} movement={:.17e} maximum_move={maximum_move:.17e} before_residual={:.17e} best_residual={best_residual:.17e} inside_movement={}", vertex.id, vertex.point, best.sub(vertex.point).length(), residual(vertex.point)?, best.sub(vertex.point).length() <= maximum_move);
let old_maximum = tolerance * 100.0 * (1.0 + vertex.point.length());
let movement = best.sub(vertex.point).length();
eprintln!("JUNCTION-OLD-MOVE vertex={} old_maximum={old_maximum:.17e} old_inside={} flip={} faces={keys:?}", vertex.id, movement <= old_maximum, (movement <= old_maximum) != (movement <= maximum_move));
}
if best.sub(vertex.point).length() <= maximum_move
&& best_residual <= residual(vertex.point)?
{
vertex.point = best;
}
}
let vertices = result
.vertices
.iter()
.map(|vertex| (vertex.id, vertex.point))
.collect::<HashMap<_, _>>();
for piece in &mut result.pieces {
let [curve_start, curve_end] = piece.curve.domain().or_refuse(KernelStage::Intersect, "domain")?;
let start = if canonical_vertices.contains(&piece.start_vertex_id) {
*vertices
.get(&piece.start_vertex_id)
.ok_or_else(|| KernelRefusal::internal(KernelStage::Intersect, "imprint.junctions", format!("missing imprint vertex {}", piece.start_vertex_id)))?
} else {
piece.curve.evaluate(curve_start).or_refuse(KernelStage::Intersect, "evaluate")?
};
let end = if canonical_vertices.contains(&piece.end_vertex_id) {
*vertices
.get(&piece.end_vertex_id)
.ok_or_else(|| KernelRefusal::internal(KernelStage::Intersect, "imprint.junctions", format!("missing imprint vertex {}", piece.end_vertex_id)))?
} else {
piece.curve.evaluate(curve_end).or_refuse(KernelStage::Intersect, "evaluate")?
};
if !canonical_vertices.contains(&piece.start_vertex_id)
&& !canonical_vertices.contains(&piece.end_vertex_id)
{
continue;
}
piece.curve = snap_curve_ends(&piece.curve, start, end)?;
for pcurve in &mut piece.pcurves {
let key = FaceKey {
operand: pcurve.operand,
face_id: pcurve.face_id,
};
let face = face_map
.get(&key)
.ok_or_else(|| KernelRefusal::internal(KernelStage::Intersect, "imprint.junctions", format!("missing imprint support face {key:?}")))?;
pcurve.pcurve = rebuild_pcurve_preserving_branch(chart_of(charts, key, face), &piece.curve, &pcurve.pcurve, emitted_pcurves.contains(&(piece.id, key)))?;
}
}
Ok(())
}
fn junction_reference(mut points: impl Iterator<Item = Vec3>) -> Vec3 {
let Some(anchor) = points.next() else {
return Vec3::default();
};
let (sum, count) = points.fold((Vec3::default(), 1usize), |(sum, count), point| {
(sum.add(point.sub(anchor)), count + 1)
});
anchor.add(sum.scale(1.0 / count as f64))
}
pub(super) fn raw_solid_extent(solid: &BrepSolid) -> f64 {
if solid.vertices.is_empty() {
return 0.0;
}
let mut low = Vec3::new(f64::INFINITY, f64::INFINITY, f64::INFINITY);
let mut high = Vec3::new(f64::NEG_INFINITY, f64::NEG_INFINITY, f64::NEG_INFINITY);
for vertex in &solid.vertices {
low.x = low.x.min(vertex.point.x);
low.y = low.y.min(vertex.point.y);
low.z = low.z.min(vertex.point.z);
high.x = high.x.max(vertex.point.x);
high.y = high.y.max(vertex.point.y);
high.z = high.z.max(vertex.point.z);
}
high.sub(low).length()
}
pub(super) fn residual_merge_bands(raw_extent: f64, tolerance: f64) -> (f64, f64) {
let band = 2.0 * assembler_weld(tolerance) * raw_extent.max(1.0);
let sep_cap = (1.5e-3 * raw_extent).max(band);
(band, sep_cap)
}
fn junction_find(parent: &mut HashMap<u64, u64>, mut x: u64) -> u64 {
loop {
let p = *parent.get(&x).unwrap_or(&x);
if p == x {
return x;
}
let grand = *parent.get(&p).unwrap_or(&p);
parent.insert(x, grand);
x = grand;
}
}
#[allow(clippy::too_many_arguments)]
pub(super) fn merge_residual_coincident_vertices(
result: &mut ImprintResultRecord,
face_map: &HashMap<FaceKey, &FaceRecord>,
edge_curves: &HashMap<(u8, u64), &EdgeRecord>,
operand_vertices: &[(u8, Vec3)],
charts: &FaceCharts<'_>,
incident_faces: &HashMap<u64, HashSet<FaceKey>>,
tolerance: f64,
raw_extent: f64,
evidence: &JunctionEvidence,
emitted_pcurves: &HashSet<(u64, FaceKey)>,
) -> Result<bool, KernelRefusal> {
let (band, sep_cap) = residual_merge_bands(raw_extent, tolerance);
let weld = assembler_weld(tolerance);
let radius_of = |id: u64| evidence.radii.get(&id).copied().unwrap_or(0.0).max(weld);
let operand_vertex_at = |point: Vec3, radius: f64| -> Option<(u8, Vec3)> {
operand_vertices
.iter()
.filter(|(_, vertex)| vertex.sub(point).length() <= radius)
.min_by(|a, b| a.1.sub(point).length().total_cmp(&b.1.sub(point).length()))
.copied()
};
let operand_guard = std::env::var("BREP_RESIDUAL_MERGE_OPERAND_GUARD").as_deref() != Ok("0");
let verts: Vec<(u64, Vec3)> = result.vertices.iter().map(|v| (v.id, v.point)).collect();
let mut parent: HashMap<u64, u64> = verts.iter().map(|(id, _)| (*id, *id)).collect();
let mut merged_point: HashMap<u64, Vec3> = verts.iter().map(|(id, p)| (*id, *p)).collect();
let mut merged_any = false;
let piece_connected: HashSet<(u64, u64)> =
if std::env::var("BREP_RESIDUAL_MERGE_PIECE_GUARD").as_deref() != Ok("0") {
result
.pieces
.iter()
.map(|piece| {
let a = piece.start_vertex_id.min(piece.end_vertex_id);
let b = piece.start_vertex_id.max(piece.end_vertex_id);
(a, b)
})
.collect()
} else {
HashSet::default()
};
for a in 0..verts.len() {
for b in (a + 1)..verts.len() {
let (id_a, p_a) = verts[a];
let (id_b, p_b) = verts[b];
let sep = p_a.sub(p_b).length();
if sep < 1e-12 || sep > sep_cap {
continue;
}
if piece_connected.contains(&(id_a.min(id_b), id_a.max(id_b))) {
continue;
}
let debug_junction = std::env::var("BREP_DEBUG_JUNCTION").is_ok();
if operand_guard {
if let (Some((operand_a, vertex_a)), Some((operand_b, vertex_b))) =
(operand_vertex_at(p_a, radius_of(id_a)), operand_vertex_at(p_b, radius_of(id_b)))
{
if operand_a == operand_b && vertex_a.sub(vertex_b).length() > weld {
if debug_junction {
eprintln!(
"residual-merge v{id_a}/v{id_b} sep={sep:.4e}: REFUSED two distinct vertices of operand {operand_a} ({:.3e} apart)",
vertex_a.sub(vertex_b).length()
);
}
continue;
}
}
}
let mut union_keys: HashSet<FaceKey> = HashSet::default();
for id in [id_a, id_b] {
if let Some(keys) = incident_faces.get(&id) {
union_keys.extend(keys.iter().copied());
}
}
if union_keys.len() < 3 {
if debug_junction {
eprintln!(
"residual-merge v{id_a}/v{id_b} sep={sep:.4e}: REFUSED union_keys={} (<3)",
union_keys.len()
);
}
continue;
}
let union_surfaces: Vec<&NurbsSurface> = union_keys
.iter()
.filter_map(|key| face_map.get(key).map(|face| &face.surface))
.collect();
if union_surfaces.len() < 3 {
continue;
}
let joint_residual = |point: Vec3| -> Result<f64, KernelRefusal> {
let mut worst = 0.0f64;
for surface in &union_surfaces {
worst = worst.max(project_point_to_surface(surface, point).or_refuse(KernelStage::Intersect, "project_point_to_surface")?.distance);
}
Ok(worst)
};
let mid = p_a.add(p_b).scale(0.5);
if joint_residual(p_a)? > band
|| joint_residual(p_b)? > band
|| joint_residual(mid)? > band
{
if debug_junction {
eprintln!(
"residual-merge v{id_a}/v{id_b} sep={sep:.4e}: REFUSED residual a={:.4e} b={:.4e} mid={:.4e} band={band:.4e} keys={}",
joint_residual(p_a)?,
joint_residual(p_b)?,
joint_residual(mid)?,
union_keys.len()
);
}
continue;
}
let mut current = mid;
let mut best = mid;
let mut best_res = joint_residual(mid)?;
for _ in 0..60 {
let mut acc = Vec3::default();
for surface in &union_surfaces {
acc = acc.add(project_point_to_surface(surface, current).or_refuse(KernelStage::Intersect, "project_point_to_surface")?.point);
}
let candidate = acc.scale(1.0 / union_surfaces.len() as f64);
let candidate_res = joint_residual(candidate)?;
if candidate_res < best_res {
best = candidate;
best_res = candidate_res;
}
if candidate.sub(current).length() <= tolerance.max(1e-12) {
break;
}
current = candidate;
}
if best_res > band
|| best.sub(p_a).length() > sep_cap
|| best.sub(p_b).length() > sep_cap
{
continue;
}
let ra = junction_find(&mut parent, id_a);
let rb = junction_find(&mut parent, id_b);
if ra == rb {
continue;
}
if std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
eprintln!(
"[merge] v{id_a} <- v{id_b} sep={sep:.3e} band={band:.3e} sep_cap={sep_cap:.3e} best_res={best_res:.3e} to=({:.7},{:.7},{:.7})",
best.x, best.y, best.z,
);
}
let rep = ra.min(rb);
let absorbed = ra.max(rb);
parent.insert(absorbed, rep);
merged_point.insert(rep, best);
merged_any = true;
}
}
if !merged_any {
return Ok(false);
}
if std::env::var("BREP_RESIDUAL_MERGE_EDGE_PLACEMENT").as_deref() == Ok("0") {
relocate_junctions(result, &mut parent, &merged_point, face_map, charts, emitted_pcurves, (4.0 * sep_cap).max(1e-3), false)?;
return Ok(true);
}
absorb_junctions_enclosed_by_merged_groups(result, face_map, incident_faces, &verts, &mut parent, band, sep_cap)?;
let retracted = place_merged_groups_on_operand_edges(result, edge_curves, operand_vertices, &verts, &mut parent, &mut merged_point, tolerance, band, evidence)?;
let mut group_size: HashMap<u64, usize> = HashMap::default();
for &(id, _) in &verts {
*group_size.entry(junction_find(&mut parent, id)).or_default() += 1;
}
let moved: HashMap<u64, Vec3> = merged_point
.iter()
.filter(|(rep, _)| group_size.get(rep).is_some_and(|size| *size > 1))
.map(|(rep, point)| (*rep, *point))
.collect();
relocate_junctions(result, &mut parent, &moved, face_map, charts, emitted_pcurves, (4.0 * sep_cap).max(1e-3), true)?;
if retracted > 0 && std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
eprintln!("[merge] retracted {retracted} edge split(s) of absorbed junctions");
}
Ok(true)
}
fn segment_on_surfaces(surfaces: &[&NurbsSurface], a: Vec3, b: Vec3, band: f64) -> Result<bool, KernelRefusal> {
for point in [a, b, a.add(b).scale(0.5)] {
for surface in surfaces {
if project_point_to_surface(surface, point).or_refuse(KernelStage::Intersect, "project_point_to_surface")?.distance > band {
return Ok(false);
}
}
}
Ok(true)
}
fn absorb_junctions_enclosed_by_merged_groups(
result: &ImprintResultRecord,
face_map: &HashMap<FaceKey, &FaceRecord>,
incident_faces: &HashMap<u64, HashSet<FaceKey>>,
verts: &[(u64, Vec3)],
parent: &mut HashMap<u64, u64>,
band: f64,
sep_cap: f64,
) -> Result<(), KernelRefusal> {
let debug = std::env::var("BREP_DEBUG_JUNCTION").is_ok();
let empty = HashSet::default();
for _ in 0..verts.len() {
let mut absorbed_any = false;
let mut members: HashMap<u64, Vec<(u64, Vec3)>> = HashMap::default();
for &(id, point) in verts {
let rep = junction_find(parent, id);
members.entry(rep).or_default().push((id, point));
}
for &(id, point) in verts {
let own_rep = junction_find(parent, id);
if members.get(&own_rep).is_some_and(|group| group.len() > 1) {
continue;
}
let mut target: Option<u64> = None;
let mut pieces = 0usize;
for piece in &result.pieces {
let other = if piece.start_vertex_id == id {
piece.end_vertex_id
} else if piece.end_vertex_id == id {
piece.start_vertex_id
} else {
continue;
};
pieces += 1;
let rep = junction_find(parent, other);
if target.is_some_and(|target| target != rep) {
target = None;
pieces = 0;
break;
}
target = Some(rep);
}
let Some(rep) = target else {
continue;
};
if pieces < 2 || rep == id {
continue;
}
let Some(group) = members.get(&rep).filter(|group| group.len() > 1) else {
continue;
};
let own_keys = incident_faces.get(&id).unwrap_or(&empty);
let mut proven: Option<(u64, f64)> = None;
for &(member_id, member_point) in group {
let sep = member_point.sub(point).length();
if sep > sep_cap {
continue;
}
let mut keys = own_keys.clone();
keys.extend(incident_faces.get(&member_id).unwrap_or(&empty).iter().copied());
let surfaces: Vec<&NurbsSurface> = keys.iter().filter_map(|key| face_map.get(key).map(|face| &face.surface)).collect();
if surfaces.is_empty() || !segment_on_surfaces(&surfaces, point, member_point, band)? {
continue;
}
if proven.is_none_or(|(_, best)| sep < best) {
proven = Some((member_id, sep));
}
}
let Some((member_id, sep)) = proven else {
continue;
};
if debug {
eprintln!("[merge] v{rep} <- v{id} enclosed: {pieces} pieces lead into the group, {sep:.3e} from member v{member_id}, band={band:.3e}");
}
parent.insert(own_rep, rep);
absorbed_any = true;
}
if !absorbed_any {
break;
}
}
Ok(())
}
fn consolidate_fuzzy_piece_ends(
result: &mut ImprintResultRecord,
face_map: &HashMap<FaceKey, &FaceRecord>,
edge_curves: &HashMap<(u8, u64), &EdgeRecord>,
charts: &FaceCharts<'_>,
tolerance: f64,
evidence: &JunctionEvidence,
emitted_pcurves: &HashSet<(u64, FaceKey)>,
) -> Result<bool, KernelRefusal> {
let weld = assembler_weld(tolerance);
let radius_of = |id: u64| evidence.radii.get(&id).copied().unwrap_or(0.0).max(weld);
let position: HashMap<u64, Vec3> = result.vertices.iter().map(|vertex| (vertex.id, vertex.point)).collect();
let mut stale: HashMap<u64, Vec3> = HashMap::default();
for piece in &result.pieces {
let [t0, t1] = piece.curve.domain().or_refuse(KernelStage::Intersect, "domain")?;
for (parameter, id) in [(t0, piece.start_vertex_id), (t1, piece.end_vertex_id)] {
let Some(point) = position.get(&id) else {
continue;
};
let gap = piece.curve.evaluate(parameter).or_refuse(KernelStage::Intersect, "evaluate")?.sub(*point).length();
if gap > tolerance && gap <= radius_of(id) {
stale.insert(id, *point);
}
}
}
if stale.is_empty() {
return Ok(false);
}
if std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
let mut ids: Vec<u64> = stale.keys().copied().collect();
ids.sort_unstable();
eprintln!("[consolidate] piece ends within their vertex's evidence radius snapped onto it: {ids:?}");
}
let mut parent: HashMap<u64, u64> = position.keys().map(|id| (*id, *id)).collect();
relocate_junctions(result, &mut parent, &stale, face_map, charts, emitted_pcurves, 0.0, true)?;
let mut retracted = 0usize;
for (id, point) in &stale {
let radius = radius_of(*id);
for record in &mut result.edge_splits {
let Some(edge) = edge_curves.get(&(record.operand, record.edge_id)) else {
continue;
};
let mut within: Vec<(f64, f64)> = Vec::new();
for ¶meter in &record.parameters {
let distance = edge.curve.evaluate(parameter).or_refuse(KernelStage::Intersect, "evaluate")?.sub(*point).length();
if distance <= radius {
within.push((distance, parameter));
}
}
if within.len() < 2 {
continue;
}
within.sort_by(|a, b| a.0.total_cmp(&b.0));
let keep = within[0].1;
let before = record.parameters.len();
record
.parameters
.retain(|¶meter| parameter == keep || !within.iter().any(|(_, other)| *other == parameter));
retracted += before - record.parameters.len();
}
}
if retracted > 0 && std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
eprintln!("[consolidate] retracted {retracted} duplicate edge split(s) within the vertex radius");
}
Ok(true)
}
fn split_edges_of(
split_points: &[((u8, u64), Vec3)],
point: Vec3,
radius: f64,
) -> HashSet<(u8, u64)> {
split_points
.iter()
.filter(|(_, split)| split.sub(point).length() <= radius)
.map(|(key, _)| *key)
.collect()
}
#[allow(clippy::too_many_arguments)]
fn place_merged_groups_on_operand_edges(
result: &mut ImprintResultRecord,
edge_curves: &HashMap<(u8, u64), &EdgeRecord>,
operand_vertices: &[(u8, Vec3)],
verts: &[(u64, Vec3)],
parent: &mut HashMap<u64, u64>,
merged_point: &mut HashMap<u64, Vec3>,
tolerance: f64,
band: f64,
evidence: &JunctionEvidence,
) -> Result<usize, KernelRefusal> {
let weld = assembler_weld(tolerance);
let radius_of = |id: u64| evidence.radii.get(&id).copied().unwrap_or(0.0).max(weld);
let mut split_points: Vec<((u8, u64), Vec3)> = Vec::new();
for record in &result.edge_splits {
let Some(edge) = edge_curves.get(&(record.operand, record.edge_id)) else {
continue;
};
for ¶meter in &record.parameters {
let point = edge.curve.evaluate(parameter).or_refuse(KernelStage::Intersect, "evaluate")?;
split_points.push(((record.operand, record.edge_id), point));
}
}
let mut groups: HashMap<u64, Vec<(u64, Vec3)>> = HashMap::default();
for &(id, point) in verts {
let rep = junction_find(parent, id);
groups.entry(rep).or_default().push((id, point));
}
let mut groups: Vec<(u64, Vec<(u64, Vec3)>)> = groups.into_iter().collect();
groups.sort_by_key(|(rep, _)| *rep);
let debug = std::env::var("BREP_DEBUG_JUNCTION").is_ok();
let mut retract: Vec<((u8, u64), Vec3, f64)> = Vec::new();
for (rep, members) in &groups {
if members.len() < 2 {
continue;
}
let mut pins: Vec<(u64, Vec3)> = Vec::new();
for &(id, point) in members {
let radius = radius_of(id);
for &(_, vertex) in operand_vertices {
if vertex.sub(point).length() <= radius && !pins.iter().any(|(_, pin)| pin.sub(vertex).length() <= weld) {
pins.push((id, vertex));
}
}
}
if pins.len() > 1 {
if debug {
eprintln!("[merge] group v{rep} ({} members) spans {} distinct operand vertices: the joint point stands", members.len(), pins.len());
}
continue;
}
let pin = pins.first().copied().filter(|&(_, vertex)| {
members.iter().all(|&(id, point)| {
point.sub(vertex).length() <= band.max(radius_of(id))
|| split_edges_of(&split_points, point, radius_of(id)).is_empty()
})
});
if pin.is_none() && !pins.is_empty() && debug {
eprintln!("[merge] group v{rep} ({} members): the operand vertex yields to a member's own edge evidence", members.len());
}
let mut holder: Option<(usize, u64, Vec3)> = pin.map(|(id, vertex)| (usize::MAX, id, vertex));
if holder.is_none() {
for &(id, point) in members {
let rank = split_edges_of(&split_points, point, radius_of(id)).len();
if rank == 0 {
continue;
}
if holder.is_none_or(|(best, best_id, _)| rank > best || (rank == best && id < best_id)) {
holder = Some((rank, id, point));
}
}
}
let Some((rank, holder_id, position)) = holder else {
continue;
};
if debug {
let joint = merged_point.get(rep).copied().unwrap_or(position);
eprintln!(
"[merge] group v{rep} ({} members) placed on v{holder_id}: {}, {:.3e} from the joint point",
members.len(),
if rank == usize::MAX { "an operand vertex".to_string() } else { format!("{rank} operand edge(s)") },
joint.sub(position).length()
);
}
merged_point.insert(*rep, position);
let holder_edges = split_edges_of(&split_points, position, radius_of(holder_id));
for &(id, point) in members {
let radius = radius_of(id);
if id == holder_id || point.sub(position).length() <= radius {
continue;
}
for key in split_edges_of(&split_points, point, radius) {
if holder_edges.contains(&key) {
retract.push((key, point, radius));
}
}
}
}
let mut retracted = 0usize;
for (key, point, radius) in retract {
let Some(edge) = edge_curves.get(&key) else {
continue;
};
let Some(record) = result
.edge_splits
.iter_mut()
.find(|record| (record.operand, record.edge_id) == key)
else {
continue;
};
let before = record.parameters.len();
record.parameters.retain(|¶meter| {
edge.curve
.evaluate(parameter)
.map(|split| split.sub(point).length() > radius)
.unwrap_or(true)
});
retracted += before - record.parameters.len();
}
Ok(retracted)
}
#[allow(clippy::too_many_arguments)]
fn relocate_junctions(
result: &mut ImprintResultRecord,
parent: &mut HashMap<u64, u64>,
merged_point: &HashMap<u64, Vec3>,
face_map: &HashMap<FaceKey, &FaceRecord>,
charts: &FaceCharts<'_>,
emitted_pcurves: &HashSet<(u64, FaceKey)>,
spur_limit: f64,
snap_moved: bool,
) -> Result<(), KernelRefusal> {
let all_ids: Vec<u64> = parent.keys().copied().collect();
let mut rep_of: HashMap<u64, u64> = HashMap::default();
for id in all_ids {
let rep = junction_find(parent, id);
rep_of.insert(id, rep);
}
let mut seen: HashSet<u64> = HashSet::default();
let mut kept = Vec::new();
for vertex in &result.vertices {
let rep = rep_of.get(&vertex.id).copied().unwrap_or(vertex.id);
if rep != vertex.id || !seen.insert(rep) {
continue;
}
let point = merged_point.get(&rep).copied().unwrap_or(vertex.point);
kept.push(ImprintVertex { id: rep, point });
}
result.vertices = kept;
let position: HashMap<u64, Vec3> =
result.vertices.iter().map(|v| (v.id, v.point)).collect();
let mut kept_pieces = Vec::new();
let mut dropped_pieces: HashSet<u64> = HashSet::default();
for mut piece in std::mem::take(&mut result.pieces) {
let start_rep = rep_of
.get(&piece.start_vertex_id)
.copied()
.unwrap_or(piece.start_vertex_id);
let end_rep = rep_of
.get(&piece.end_vertex_id)
.copied()
.unwrap_or(piece.end_vertex_id);
let repointed = start_rep != piece.start_vertex_id || end_rep != piece.end_vertex_id;
let moved = repointed
|| (snap_moved && (merged_point.contains_key(&start_rep) || merged_point.contains_key(&end_rep)));
piece.start_vertex_id = start_rep;
piece.end_vertex_id = end_rep;
if moved {
let start = *position
.get(&start_rep)
.ok_or_else(|| KernelRefusal::internal(KernelStage::Intersect, "imprint.junctions", format!("missing merged imprint vertex {start_rep}")))?;
let end = *position
.get(&end_rep)
.ok_or_else(|| KernelRefusal::internal(KernelStage::Intersect, "imprint.junctions", format!("missing merged imprint vertex {end_rep}")))?;
if repointed && start_rep == end_rep && curve_length_rough(&piece.curve)? < spur_limit {
dropped_pieces.insert(piece.id);
continue;
}
let [t0, t1] = piece.curve.domain().or_refuse(KernelStage::Intersect, "domain")?;
let start_gap = piece.curve.evaluate(t0).or_refuse(KernelStage::Intersect, "evaluate")?.sub(start).length();
let end_gap = piece.curve.evaluate(t1).or_refuse(KernelStage::Intersect, "evaluate")?.sub(end).length();
if !repointed && start_gap <= 1e-12 && end_gap <= 1e-12 {
kept_pieces.push(piece);
continue;
}
piece.curve = snap_curve_ends(&piece.curve, start, end)?;
if std::env::var("BREP_DEBUG_JUNCTION").is_ok() {
let [t0, t1] = piece.curve.domain().or_refuse(KernelStage::Intersect, "domain")?;
let left_start = piece.curve.evaluate(t0).or_refuse(KernelStage::Intersect, "evaluate")?.sub(start).length();
let left_end = piece.curve.evaluate(t1).or_refuse(KernelStage::Intersect, "evaluate")?.sub(end).length();
eprintln!(
"[relocate] piece {} v{start_rep}->v{end_rep} gaps {start_gap:.3e}/{end_gap:.3e} -> {left_start:.3e}/{left_end:.3e}{}",
piece.id,
if left_start > 1e-9 || left_end > 1e-9 { " NOT SNAPPED (beyond the snap cap)" } else { "" }
);
}
for pcurve in &mut piece.pcurves {
let key = FaceKey {
operand: pcurve.operand,
face_id: pcurve.face_id,
};
let face = face_map
.get(&key)
.ok_or_else(|| KernelRefusal::internal(KernelStage::Intersect, "imprint.junctions", format!("missing imprint support face {key:?}")))?;
pcurve.pcurve = rebuild_pcurve_preserving_branch(
chart_of(charts, key, face),
&piece.curve,
&pcurve.pcurve,
emitted_pcurves.contains(&(piece.id, key)),
)?;
}
}
kept_pieces.push(piece);
}
result.pieces = kept_pieces;
if !dropped_pieces.is_empty() {
for face in &mut result.by_face {
face.piece_ids.retain(|id| !dropped_pieces.contains(id));
}
}
Ok(())
}