1use crate::surface_surface_intersection::SurfaceIntersectionCurve;
51use crate::tangent_curve::classify_tangent_contact_at;
52use crate::{NurbsSurface, Vec3};
53
54const NODE_WITNESS_CROSS: f64 = 0.25;
61
62pub fn crossing_node_split_enabled() -> bool {
65 crate::switch!("BREP_SSI_NODE_SPLIT", true)
66}
67
68#[derive(Clone, Copy, Debug)]
71struct Node {
72 point: Vec3,
73 params_a: [f64; 2],
74 params_b: [f64; 2],
75 radius: f64,
76}
77
78fn periods(surface: &NurbsSurface) -> Result<[f64; 2], String> {
80 let (closed_u, closed_v) = surface.closed_directions()?;
81 let [u0, u1] = surface.domain_u()?;
82 let [v0, v1] = surface.domain_v()?;
83 Ok([
84 if closed_u { u1 - u0 } else { 0.0 },
85 if closed_v { v1 - v0 } else { 0.0 },
86 ])
87}
88
89fn continue_params(params: [f64; 2], near: [f64; 2], periods: [f64; 2]) -> [f64; 2] {
91 let mut out = params;
92 for index in 0..2 {
93 let period = periods[index];
94 if period > 0.0 {
95 out[index] += period * ((near[index] - params[index]) / period).round();
96 }
97 }
98 out
99}
100
101fn crossing(first: &NurbsSurface, second: &NurbsSurface, a: [f64; 2], b: [f64; 2]) -> Option<f64> {
102 let normal_a = first.normal(a[0], a[1]).ok()?;
103 let normal_b = second.normal(b[0], b[1]).ok()?;
104 Some(normal_a.cross(normal_b).length())
105}
106
107fn find_nodes(
109 first: &NurbsSurface,
110 second: &NurbsSurface,
111 curves: &[SurfaceIntersectionCurve],
112 tolerance: f64,
113) -> Result<Vec<Node>, String> {
114 let mut nodes: Vec<Node> = Vec::new();
115 for curve in curves {
116 let count = curve.points.len();
117 if count < 3 || curve.tangent || curve.params_a.len() != count || curve.params_b.len() != count {
118 continue;
119 }
120 let crossings: Vec<f64> = (0..count)
121 .map(|index| {
122 crossing(first, second, curve.params_a[index], curve.params_b[index]).unwrap_or(f64::INFINITY)
123 })
124 .collect();
125 for index in 0..count {
126 let value = crossings[index];
127 if !(value < NODE_WITNESS_CROSS) {
128 continue;
129 }
130 let before = if index > 0 { crossings[index - 1] } else if curve.closed { crossings[count - 2] } else { f64::INFINITY };
131 let after = if index + 1 < count { crossings[index + 1] } else if curve.closed { crossings[1] } else { f64::INFINITY };
132 if value > before || value > after {
133 continue;
134 }
135 let witness = curve.points[index];
136 if nodes.iter().any(|node| node.point.sub(witness).length() <= node.radius) {
137 continue;
138 }
139 let Some(contact) = classify_tangent_contact_at(first, second, witness, tolerance)? else {
140 continue;
141 };
142 if !contact.saddle || contact.is_extended() || contact.curvature_scale <= 0.0 {
143 continue;
144 }
145 let before = if index > 0 { index - 1 } else if curve.closed { count - 2 } else { index };
148 let after = if index + 1 < count { index + 1 } else if curve.closed { 1 } else { index };
149 let chord = curve.points[before]
150 .sub(witness)
151 .length()
152 .max(curve.points[after].sub(witness).length());
153 let reach = (1e-2 * (1.0 + witness.length())).max(chord);
154 if contact.point.sub(witness).length() > reach {
155 continue;
156 }
157 if nodes.iter().any(|node| node.point.sub(contact.point).length() <= tolerance * 10.0) {
158 continue;
159 }
160 nodes.push(Node {
161 point: contact.point,
162 params_a: contact.params_a,
163 params_b: contact.params_b,
164 radius: 3.0 * (2.0 * tolerance / contact.curvature_scale).sqrt(),
165 });
166 }
167 }
168 Ok(nodes)
169}
170
171#[derive(Clone, Copy, Debug)]
175struct Cut {
176 left: usize,
177 right: usize,
178 node: usize,
179}
180
181fn point_segment_distance(point: Vec3, start: Vec3, end: Vec3) -> f64 {
182 point.distance_to_segment(start, end, crate::tolerance::FP_FLOOR)
183}
184
185fn passages(points: &[Vec3], cyclic: bool, nodes: &[Node]) -> Vec<Cut> {
195 let count = points.len() as isize;
196 let at = |index: isize| -> Vec3 { points[index.rem_euclid(count) as usize] };
197 let mut cuts: Vec<(isize, isize, usize)> = Vec::new();
199 for (node_index, node) in nodes.iter().enumerate() {
200 let distance = |index: isize| at(index).sub(node.point).length();
201 for index in 0..count {
202 if !cyclic && (index == 0 || index + 1 == count) {
203 continue;
206 }
207 let here = distance(index);
208 if here > distance(index - 1) || here >= distance(index + 1) {
209 continue;
210 }
211 let before = at(index - 1);
212 let after = at(index + 1);
213 let near_chord = point_segment_distance(node.point, before, at(index))
214 .min(point_segment_distance(node.point, at(index), after));
215 let longest = before.sub(at(index)).length().max(after.sub(at(index)).length());
216 if near_chord > node.radius.max(0.25 * longest) {
217 continue;
218 }
219 let (left, right) = if here <= node.radius {
220 let mut low = index;
223 let mut high = index;
224 let limit = count / 2 - 1;
225 while index - low < limit && distance(low - 1) <= node.radius && (cyclic || low > 1) {
226 low -= 1;
227 }
228 while high - index < limit && distance(high + 1) <= node.radius && (cyclic || high + 2 < count) {
229 high += 1;
230 }
231 (low - 1, high + 1)
232 } else if at(index).sub(node.point).dot(after.sub(before)) < 0.0 {
233 (index, index + 1)
234 } else {
235 (index - 1, index)
236 };
237 if !cyclic && (left < 0 || right >= count) {
238 continue;
239 }
240 cuts.push((left, right, node_index));
241 }
242 }
243 cuts.sort_by_key(|cut| cut.0);
244 let mut merged: Vec<(isize, isize, usize)> = Vec::new();
246 for cut in cuts {
247 if let Some(last) = merged.last() {
248 if cut.0 < last.1 {
249 continue;
250 }
251 }
252 merged.push(cut);
253 }
254 if cyclic && merged.len() > 1 {
255 let first = merged[0];
256 let last = merged[merged.len() - 1];
257 if last.1 > first.0 + count {
258 merged.pop();
259 }
260 }
261 merged
262 .into_iter()
263 .map(|(left, right, node)| Cut {
264 left: left.rem_euclid(count) as usize,
265 right: right.rem_euclid(count) as usize,
266 node,
267 })
268 .collect()
269}
270
271#[allow(clippy::too_many_arguments)]
274fn piece(
275 curve: &SurfaceIntersectionCurve,
276 indices: &[usize],
277 start_node: Option<&Node>,
278 end_node: Option<&Node>,
279 periods_a: [f64; 2],
280 periods_b: [f64; 2],
281) -> Option<SurfaceIntersectionCurve> {
282 if indices.is_empty() {
283 return None;
284 }
285 let mut points: Vec<Vec3> = indices.iter().map(|&index| curve.points[index]).collect();
286 let mut params_a: Vec<[f64; 2]> = indices.iter().map(|&index| curve.params_a[index]).collect();
287 let mut params_b: Vec<[f64; 2]> = indices.iter().map(|&index| curve.params_b[index]).collect();
288 let mut tangents: Option<Vec<Vec3>> =
289 curve.tangents.as_ref().map(|tangents| indices.iter().map(|&index| tangents[index]).collect());
290 if let Some(node) = start_node {
291 let a = continue_params(node.params_a, params_a[0], periods_a);
292 let b = continue_params(node.params_b, params_b[0], periods_b);
293 let chord = points[0].sub(node.point).normalized().ok();
294 points.insert(0, node.point);
295 params_a.insert(0, a);
296 params_b.insert(0, b);
297 if let Some(tangents) = tangents.as_mut() {
298 let first = tangents[0];
299 tangents.insert(0, chord.unwrap_or(first));
300 }
301 }
302 if let Some(node) = end_node {
303 let last = points.len() - 1;
304 let a = continue_params(node.params_a, params_a[last], periods_a);
305 let b = continue_params(node.params_b, params_b[last], periods_b);
306 let chord = node.point.sub(points[last]).normalized().ok();
307 points.push(node.point);
308 params_a.push(a);
309 params_b.push(b);
310 if let Some(tangents) = tangents.as_mut() {
311 let previous = tangents[last];
312 tangents.push(chord.unwrap_or(previous));
313 }
314 }
315 if points.len() < 2 {
316 return None;
317 }
318 Some(SurfaceIntersectionCurve {
319 points,
320 params_a,
321 params_b,
322 closed: false,
323 tangents,
324 span_deviation: curve.span_deviation,
325 tangent: false,
326 })
327}
328
329pub fn split_at_crossing_nodes(
333 first: &NurbsSurface,
334 second: &NurbsSurface,
335 curves: Vec<SurfaceIntersectionCurve>,
336 tolerance: f64,
337) -> Result<(Vec<SurfaceIntersectionCurve>, Vec<Vec3>), String> {
338 let nodes = find_nodes(first, second, &curves, tolerance)?;
339 if nodes.is_empty() {
340 return Ok((curves, Vec::new()));
341 }
342 let periods_a = periods(first)?;
343 let periods_b = periods(second)?;
344 let debug = crate::debug_flag!("BREP_DEBUG_NODE_SPLIT");
345 let mut used = vec![false; nodes.len()];
349 let mut out = Vec::with_capacity(curves.len() + 2 * nodes.len());
350 for curve in curves {
351 let full = curve.points.len();
352 if full < 3 || curve.tangent || curve.params_a.len() != full || curve.params_b.len() != full {
353 out.push(curve);
354 continue;
355 }
356 let cyclic = curve.closed && curve.points[0].sub(curve.points[full - 1]).length() <= tolerance * 10.0;
358 let count = if cyclic { full - 1 } else { full };
359 let cuts = passages(&curve.points[..count], cyclic, &nodes);
360 if debug {
361 eprintln!(
362 "node split: branch of {full} stations (closed {cyclic}), nodes {:?}, cuts {:?}",
363 nodes.iter().map(|node| (node.point, node.radius)).collect::<Vec<_>>(),
364 cuts
365 );
366 }
367 if cuts.is_empty() {
368 out.push(curve);
369 continue;
370 }
371 for cut in &cuts {
372 used[cut.node] = true;
373 }
374 let mut pieces = Vec::new();
375 if cyclic {
376 for (index, cut) in cuts.iter().enumerate() {
377 let next = cuts[(index + 1) % cuts.len()];
378 let mut indices = Vec::new();
379 let mut station = cut.right;
380 let stop = (next.left + 1) % count;
381 loop {
382 indices.push(station);
383 station = (station + 1) % count;
384 if station == stop || indices.len() >= count {
385 break;
386 }
387 }
388 pieces.push(piece(&curve, &indices, Some(&nodes[cut.node]), Some(&nodes[next.node]), periods_a, periods_b));
389 }
390 } else {
391 let mut start = 0usize;
392 let mut start_node: Option<&Node> = None;
393 for cut in &cuts {
394 let indices: Vec<usize> = (start..=cut.left).collect();
395 pieces.push(piece(&curve, &indices, start_node, Some(&nodes[cut.node]), periods_a, periods_b));
396 start = cut.right;
397 start_node = Some(&nodes[cut.node]);
398 }
399 let indices: Vec<usize> = (start..count).collect();
400 pieces.push(piece(&curve, &indices, start_node, None, periods_a, periods_b));
401 }
402 let pieces: Vec<SurfaceIntersectionCurve> = pieces.into_iter().flatten().collect();
403 if debug {
404 eprintln!(
405 "node split: -> {} piece(s) of {:?} stations",
406 pieces.len(),
407 pieces.iter().map(|piece| piece.points.len()).collect::<Vec<_>>()
408 );
409 }
410 out.extend(pieces);
411 }
412 Ok((
413 out,
414 nodes
415 .iter()
416 .zip(&used)
417 .filter(|(_, used)| **used)
418 .map(|(node, _)| node.point)
419 .collect(),
420 ))
421}
422